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ABSTRACT 

We show that observations of high-redshift Lya emitters (LAEs) have the poten- 
tial to provide definitive evidence for reionization in the near future. Using 200 Mpc 
radiative transfer simulations, we calculate the effect that patchy reionization has on 
the line profile, on the luminosity function, and, most interestingly, on the cluster- 
ing of emitters for several realistic models of reionization. Reionization increases the 
measured clustering of emitters, and we show that this enhancement would be essen- 
tially impossible to attribute to anything other than reionization. Our results motivate 
looking for the signature of reionization in existing LAE data. We find that for stellar 
reionization scenarios the angular correlation function of the 58 LAEs in the Subaru 
Deep Field z = 6.6 photometric sample is more consistent with a fully ionized uni- 
verse (mean volume ionized fraction ~ 1) than a universe with Xi < 0.5 at > 2-a 
confidence level. Measurements in the next year on Subaru will increase their z = 6.6 
LAE sample by a factor of five and tighten these limits. If the clustering signature of 
reionization is detected in a LAE survey, a comparison with a Lyman-break or a Ha 
survey in the same field would confirm the reionization hypothesis. We discuss the op- 
timal LAE survey specifications for detecting reionization, with reference to upcoming 
programs. 

Key words: cosmology: theory - diffuse radiation - intergalactic medium - large- 
scale structure of universe - galaxies: high redshift - line: profiles 



1 INTRODUCTION 

The reionization epoch - when the hydrogen in the Universe 
was ionized by photons produced in the first galaxies - re- 
, mains one of the least explored periods in cosmology. There 
" is substantial uncertainty as to when reionization occurred. 
The latest cosmic microwave background (CMB) data show 
that the mean redshift of reionization was Zrei — 1 1 zb o , 
but are consistent with Zrei = at 3-cr (|Page et all 120061 ). 
and measurements of Lya absorption in the spectra of high- 
redshift quasars sug gest that reionization ended at z ~ 6 
(|Becker et all 1200 II : I White et al.1 12003: Fan et al.|[2006l ). 

Whether observations have already pinned down the 
redshift of bubble percolation to be z ~ 6 using the Lya 
absorption features in the spectra of high-redshift quasars - 
the Lya forest - is under debate. The appearance of Gunn- 
Peterson absorption troughs in the z > 6 Lya forest may 
signify a ch ange in the ionization state of the i ntergalac- 
tic m edium (jBecker et al.ll200ll : lwhite et al.|[2003l : lFan et all 
I2OO6I I. In addition, the rapid decrease of Lya forest trans- 
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mission at z > 6 (as well as the decrease in transmission 
in the Ly/3 and the Ly7 forests), the large variance in the 
transmission between sight-lines, and the sizes of the prox- 
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However, owing to the resonant nature of the Lya line, 
Lya transmission is essentially zero if the gas is more neu- 
tral than one part in a thousand, which makes the Lya for- 
est insensitive to the order unity fluctuations in the ion- 
ized fraction that define reionization. Because of this limi- 
tation, arguments for t he end of reionizatio n based on the 
forest are controversial. iBecker et al.l (|2006l l proposed that 
the rapid evolution of Lya transmission in the z ~ 6 forest 
is consisten t with an ionized in t ergala ctic medium (IGM). 
iLidz et al.1 (|2006h and lLiu et all (|2006l l showed that the ob- 
served variance in the transmission between sight-lines is 
consistent with the variance one ex pect s from an i onized 
IGM, and Bolton k Haehnelt (2007) and' Lidz etHI (|200/t ) 
demonstrated that current data cannot distinguish the prox- 
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imity region owing to enhanced ionizing flux from the effect 
of an HII region. 

High-redshift QSOs will probably not reveal additional 
information about reionization in the next few years, as all- 
sky surveys that observe in the near-infrared are required 
to push to higher redshifts. The Lya forest in the spectrum 
of high-redshift gamma ray bursts (GRBs) holds more im- 
mediate promise, but the presence of damped Lya systems 
within GRB host galaxies complicates the determination 
of the ionization state of the IGM. Despite this difficulty, 
through measurements of t he Lyg and th e Ly/3 forests in the 
spectrum of a single GRB, Totani et al. (2006") derived the 
constraint Xi > 0.4 at z = 6.3, assuming a uniformly ionized 
IGM, where Xi is the global ionized fraction of the IGM. Ac- 
counting for a patchy reionization process will weaken this 
limit. 

Observations of high-redshift 21cm emission have the 
potential to put the strongest constraints on Xi in the long 
term. In principle, such studies can image the neutral hy- 
drogen in the Universe as it became ionized. However, 21cm 
observations must first overcome terrestrial radio broadcasts 
as well as galactic foregrounds that are four orders of mag- 



priate observing strategy - will effectively be able to image 
HII regions during reionization. 

There are several surveys that have targeted or 
are targeting LAEs at epochs whe n the Universe may 



have been significantly neutral (e.g., 'Kodaira et al.' '2003 



nitude larger than the signal (e.g., [Zal darriaga et al. I l2004l : 



iMcQuinn et al.ll2006l : fFurlanetto et al.l l2006). Despite these 
difficulties, the Mileura Widefield Array and the Low Fre- 
quency Array will start de dicated observations of this emis- 
sion in the coming years ([Morales k HewittI l2004l : Ide VosI 
l2004l ). and the 21cm Array and the Giant Meter- wave Radio 
Telesco pe are already b eing used to perform these measure- 
ments (" Pen et al.ir200 5V 

Studies of the CMB will never be able to provide a 
direct detection of reionization by imaging individual HII 
regions. However, measurements of CMB anisotropics can 
improve constraints on the global properties of reioniza- 
tion. A more precise measurement of CMB polarization 
anisotropics at low multipoles has the potential to pro- 
vide some information about the duration of reionization 
in addition to tighter constraints on the average redshift 
of reionization (e.g., Keating & Miller 2006). The Planck 
mission will measure the large-scale E-mode polarization 
anisotropics to nearly the cosmic- variance limit in the next 
few years. In addition, a sizable fraction of arc-minute-scale 
CMB anisot ropics are imprinted during patchy reionization 
(|Santos et al . 2003; Zahn et al. 2005; McQuinn et al. 2005), 
and the Atacama Cosmology Telescope and the South Pole 
Telescope are starting to observe these anisotropics Q These 
observations have the potential to provide information about 
the duration of reionization and about the sizes of the HII 
regions. 

In this paper, we demonstrate that Lya emitter (LAE) 
surveys will likely offer the first irrefutable proof as to pre- 
cisely when reionization occurred, at least if the Universe is 
significantly neutral for z < 7 as the high-redshift quasars 
suggest. We show that the z = 6.6 LAE sample from Sub- 
aru can already put Umits on Xi with just 58 photometrically 
confirmed LAEs. Surveys on existing telescopes will also put 
constraints on the neutral fraction at higher redshifts, and 
the James Webb Space Telescope ( JWST) - with the appro- 
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iRhoad s et al. 2004; Kashikawa et al. 2006; Cubv et al. 20071 : 
Willis & Courbin 2005; Sta rk et al.ll2007|). and other pro- 
grams will begin takin g data soon (e .g.JCasah et al.ll2006l : 
iMcPherson et al.1 l2006l : iHorton et all l2004l l. Thesesurvevs 
take advantage of the fact that a gala xy can produce a siz- 
able fraction of its flux in the Lya Une ((Partridge PeeblesI 
1967). The transmission of the Lya line from a galaxy is de- 
creased if the neighborhood of an emitter is largely neutral, 
and this modulation can be used to probe the epoch of reion- 
ization ("Miralda-Escude'igQS"; Madau k Rees 2000;'Haiman| 
2002; Santos 2004; Furlanetto et al. 2004c; Furlanetto et al. 
l2006i : Malhotra k Rhoadsl 12006a). 

Most previous theoretical studies of high-redshift LAEs 
have focused on using the evolution in the LAE luminosity 
function or in the average Lya line profile to constrain reion- 
ization. The first predictions for the impact of reionization 
on the LAE luminosity function and line profile were based 
on the assumption that each emitter sits in its own HII re- 
gion until the final bubble percolation stage, at which time 
all of th e HII regions quickly mer ge and reionization is com- 
pleted (jSantod Eooi : iHaiman "k G en 2005). In reality, the 
HII regions can be much larger throughout reionization than 
this simple model pr edicts owing to the c l ustering of high- 
redshi ft sources (e.g. JSokasian et "al]|2003l : iFurlanetto et al.l 
I2004al li3). When the Universe has Xi = 0.5, typically a LAE 
will be in an HII region created by thousands of galax- 
ies. Because of this clustering, both the evolution of the 
LAE luminosity function and the change in the shape of 
the average line profile of an emitter are a weaker func- 
tion of the neutral fraction than predicte d by models based 
on a single source per HII region (.Furlanetto et al. I l2004d : 
IFurlanetto et al ] |2006l : lMcQuinn et al.ll20oi L 

The measured z = 6.5 LAE luminosity function has al- 
ready been used to constrain Xi. Malhotra k Rhoads (20^ 
argued that the data is consistent with no evolution in the 
luminosity function bet ween z — 5.7 and z = 6.5, a nd this 
lack of evolution allowed iMalhotra k Rhoads! (|2006b[ ) to de- 
rive the limit > 0.2 at z = 6 .5. However, in th e newest and 
largest LAE survey at z > 6, iKashikawa et al.l (I2OO6) finds 
significant evolution between the z = 5.7 and the z = 6.6 
luminosity functions for luminosities > 2 x 10^^ erg s~^, and 
they suggest that thi s evolution could be evidence for reion- 
ization. Furthermore, 'lye et al." (2007') finds that the number 
density of LAEs with luminosities above 1 x 10^^ erg s~^ at 
2; = 7 is 20% — 40% smaller than this number density at 
z = 6.6. 

In addition to the luminosity function and the Lya line 
profile, the HII bubbles will influence the measured cluster - 
ing of LAEs (Furlanetto et al.ll2006l : iMcQuinn eral]l2006l ). 
Detecting reionization through its impact on observed clus- 
tering will be the most fool-proof method to identify pockets 
of neutral hydrogen in the IGM with LAEs, since astro- 
physical uncertainties and observational systematics cannot 
induce large-scale correlations similar to ^ 10 Mpc HII re- 
gions during reionization. We show that a LAE survey using 
clustering measurements can confirm whether reionization 
is occurring at z = 6.6 in the immediate future and that 
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upcoming surveys will place constraints on reionization at 
higher redshifts. 

This paper is organized as follows. In Section [21 we de- 
scribe the simulations of reionization that we employ as well 
as the post-processing calculations used to generate mock 
Lya surveys from the simulation outputs. These mock LAE 
surveys are featured in Section [3] The effect of reionization 
on the luminosity function is discussed in Section |4] (and 
on the line profile in Appendix |A}. However, the primary 
focus of this work is to investigate the impact that reioniza- 
tion has on the measured clustering properties of LAEs and 
to quantify the detectability of reionization-induced cluster- 
ing in present surveys (Section [5] and, for void statistics. 
Appendix [BJ . We conclude in Section [6] with a discussion 
of the capabilities of upcoming LAE surveys to constrain 
reionization. 



2 CALCULATIONS 
2.1 Simulations 

We use two 1024^ N-body simulations generated with the 
TreePM code L-Gadget-2 (^ringel 2005) to model the den- 
sity field, one in a box of size 94Mpc with outputs every 
50 Myr and the other in a box of size 186 Mpc with out- 
puts every 20 Myr. A Friends-of-Friends algorithm with a 
linking length of 0.2 times the mean inter-particle spacing 
is used to identify halos. We employ the higher-resolution 
94 Mpc box to study Une properties. This box is run us- 
ing a ACDM cosmology with = 1, as = 0.9, — 0.3, 

— 0.7, Qb — 0.04, a nd h — 0.7, consist ent with the 
first year WMAP results (Spergel et al.l l2003h . We use the 
186 Mpc box to study the luminosity function and clustering 
properties of emitters. This volume provides a larger sam- 
ple of the structures that are present during reionization 
than does the smaller box. The 186 Mpc box is run with 
a ACDM cosmology updated to be m ore consistent with 
the latest WMAP results (Spe rgel et al. 2006) with = 1, 
(78 0.8, 0.27, Qa 0.73, 0.046, and h 0.7. 

The assumptions about the cosmology do not significantly 
aff ect the morphology of reionization if we compare at fixed 
Xi (|McQ uinn et al. 2006). 

The halo mass f unction measured f rom t he N-body sim- 
ulations matches the lSheth TormenI (|2002l l mass function 
for groups with 32 or more particles. Thirty- two particle 
groups correspond to a halo mass of m = 1 X 10^ Mo m 
the 94 Mpc box and m = 8 x 10^ M© in the 186 Mpc box. 
We would like to resolve halos down to the atomic hydrogen 
cooling mass mcooi - halos with virial temperature of 10^ K 
- which corresponds to the minimum mass galaxy that can 
form stars (mcooi ~ 10^ at ^ = 8). We add unresolved 
halos (halos that would be comprised of fewer than 32 par- 
ticles) into the simulation using a Press- Schechter merger 
tree. This algorithm is similar to the PThalo method for gen - 
erating mock catalogs of galaxies (|Sheth Lemsonlll999[ ), 
and accounts for both density and Poisson fluctuations in 
the a bundance of halos. We demonstrate in lMcQuinn et aTl 
([20061 1 that this method reproduces well the power spectrum 
and mass function of dark matter hal os seen in simulations . 

We use an improved version of the Sokasian et al ^ ("2001', 
I2OO3I . HqqI) and iMcQuinn et~all (j2006) radiative transfer 



code. We have refined the algorithm to group sources more 
efficiently for our work, and this improvement is discussed in 
Appendix [D] Our approach is optimized to simulate reion- 
ization, making several justified simplifications to speed up 
the computation. The code takes the gridded density field 
(generated from the N-body simulation assuming that the 
baryons trace the dark matter) and a list of the ionizing 
sources as input, and it casts rays from each source to com- 
pute the ionization field. We assume that the sources have 
a soft UV spectrum that scales as z/~^ (consistent with a 
POPII IMF), which is used to calculate the photo- ionization 
state of the gas. The radiative transfer code assumes per- 
fectly sharp HIT fronts, tracking the front position at sub- 
grid scaleso Sharp fronts are an excellent approximation for 
sources with a soft UV spectrum, where the mean free path 
for ionizing photons is kiloparsecs, substantially smaller than 
the cell size in the radiative transfer calculations. 

The radiative transfer calculations presented in this pa- 
per were performed on a 256^ grid for the 94 Mpc box and 
on a 512*^ grid for the 186 Mpc box. The radiative transfer 
calculations took between three days and two weeks on a 2.2 
GHz AMD Opteron processor. 

2.2 Lya Line Calculation 

The shape and energetics of a Lya line from a galaxy are 
determined by complex processes, many of which cannot 
be modeled accurately. The Lya emission depends on the 
amount of dust in the galaxy, the velocity profile of the gas 
in the halo, and the fraction of ionizing photons that escape 
from the galaxy and influence the photo-ionization sta te of 
the nea rby IGM (Hansen & Oh 2006; Diikstra et al. 20071 : 
ISantos et al.1 l2003l : ITapken et al.1 120071 : ITasitsiomil l26o6h . 
Statistics that can be used to isolate the effect that reion- 
ization has on the emitters from all of the complicated, un- 
certain astrophysics are essential to probe reionization with 
LAEs. In order to construct mock surveys from which we 
can measure various statistics, we must make simplifying 
assumptions about the Lya emission that escapes from the 
vicinity of a galaxy. 

To calculate the emitter line profile and its transmission 
we employ two schemes, a detailed, expensive method and 
a fast, relatively inexpensive technique. 

Method 1: We assume a Gaussian intrinsic line profile - 
the line profile that escapes from the vicinity of the galaxy - 
with its width set by the circular velocity at the virial radius 
of the emitter host halo (at, — 1/^ v-vir/c). This intrinsic line 
profile is wh at is commonly assumed in the literature (e.g., 
ISantodliooi and lPiikstra et al.|[2007l ). but the outgoing line 
can be significantly more complicated ([Ha nsen & Oh 20061 : 
Diiks tra et aD I2OO6I : iTapken et all l2007i r[Tasitsiomil l2006l ). 
See the discussion in Appendix|X]for the physical motivation 
for the assumed line profile. In the absence of a galactic wind, 
we assume the intrinsic line is centered at Lya. Given the 
intrinsic line from an emitter, we compute the absorption 
along a skewer through the simulation box, using the density, 
velocity, and photo-ionization state of the gas. The radiative 
transfer simulations do not compute the temperature of the 

This is not true for self-shielded regions, which can remain neu- 
tral behind the front. 
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gas self-consistently. For simplicity, we take the gas in all 
regions to be at lO^iT. The value of the gas temperature 
and the amplitude of the sub-grid-scale density fluctuations 
do not influence the amount of damping wing absorption - 
the effect that is most relevant to this study. 

Method 2: First, we set the residual neutral fraction 
within each HII region to zero. The residual neutral gas pri- 
marily influences frequencies blueward of Lya, which this 
method assumes are fully absorbed. Next, we calculate the 
damping wing optical depth to absorption along a ray 
through the IGM, ignoring the effect of pecuUar velocities. 
We do this for only a single frequency that starts off at 
y' — (1 — / c) for each source, where is the line cen- 
ter and is the wind velocity (which is typically taken to be 
zero)If| We then set Lobs, a — %x La exp[— Ta(z^O], where 7^ 
accounts for resonant absorption and is typically between 
0.1 and 0.5 in the absence of outflows (|Santos et al.l l200"3l 
and Appendix [A|. The exact value of 7^ or its scaling with 
halo mass is not important for our present study since 
is degenerate with the luminosity of the emitter prior to 
absorption in the IGM, La. This method is in contrast to 
Method 1 in which for each source we calculate Ta at many 
different frequencies. 

It turns out that both methods for computing the line 
profile result in very similar results. (See Figure IA2I for a 
comparison of the two methods.) The agreement between 
the two methods occurs because for observed emitters the 
virial velocity of the halo, which sets the intrinsic width of 
the line in our model, is typically much smaller than the 
relative velocity between the LAE and the HII front. The 
agreement between the two methods justifies the use of the 
cheaper method. Method 2, for computations in which the 
shape of the line profile is not required. 

For reference, the radius of a bubble that produces a 
damping wing absorption cross section of unity at source- 
frame frequency u ^ uq is (|Loeb Rvbickilll999 ) 



Rb ~ 1.1 Xf 



0.046 a 



pMpc + 



1^0 



1^0 H{Zg 



(1) 



where pMpc denotes proper Mpc, xh is the neutral fraction 
outside of the bubble, and Zg is the redshift of a LAE in the 
center of the bubble. Equation ([1]) is correct in detail only 
for a homogeneous xh. In practice, the IGM is far from 
uniform and the effective xh along each line of sight is not 
equal to the globally averaged xh- Also note that Tai^^) ~ 
900 kms-' XH [(1 + %)/8]^/^ [H{zg)Rb - c{v - - 
the optical depth scales i nversely with the bubble radius. 
(See 'Mirald a-Escudel ll998' for the exact expression for Ta as 
a function of the bubble size.) If a LAE produces its own 
HII bubble then the size of this bubble is 



Rh = 0.15 



S 



tage fesc 

Mq yr-1 108 yr OT 



1 + z 



pMpc, (2) 



where S is the star formation rate, which is estimated 
to be ap proxim ately 1-10 Mq yr~^ for observed emitters 
(|Taniguchi et al.ll2005l ). and tage is the lifetime of the emit- 
ter. Equation (|2]) does not account for recombinations and 



^ Note that the actual wind velocity is probably closer to 0.5 -Uw 
such that v^Jv characterizes the shift of the Lyo; line rather than 
the wind velocity. 



assumes a Salpeter IMF with 1.5 x 10^^ 5'/(M0yr"^) ion- 
izing photons s" (|Hui et all liooi ). Note that typically 
Ta(z^o) > 1 for an emitter that sits in an HII region cre- 
ated entirely by itself. In CDM theory, high redshift galaxies 
are very clustered such that the situation of a single source 
creating an HII region is a rarity. 



2.3 Simple Model for Intrinsic Distribution of 
LAEs 

Given the methods discussed in Section [2?2l for calculating 
the line profile, we also need an algorithm that tells us where 
the LAEs are located to generate mock surveys. For the 
analysis in this paper, we adopt a simple model to populate 
the halos in the simulation with LAEs. We describe a mock 
survey of LAEs with three parameters: 

• ^min - minimum mass halo that hosts an observed LAE 

• /e - fraction of halos that host LAEs 

• Fc - fraction of the LAE sample contaminated by lower- 
redshift interlopers 

We assume that each halo in the simulation box has at most 
one LAE, which is reasonable at z > 6 since the gas cooling 
time in all halos is shorter than a Hubble time. We also 
assume that /e is not a function of h alo mass|3 This model 
is simil ar to those for LA Es use d in [ Piikstra et al. J 2006r ) 
and in IStark et al.1 (|2007f ). Both IPiikstra et al.1 (|2006h and 
IStark et al.l (|200'^ find that this simple parameterization is 
able to provide adequate fits to the z = 5.7 and z = 6.6 
luminosity functions even when fixing Fe- 
in the context of this model, we can write expressions 
for the intrinsic number density of observed emitter can- 
didates nE = /e (1 + Fc)nh(rnmin) as well as the intrin- 
sic 3-D n£;-weighted power spectrum of emitter candidates 
Pint,E(k) = (l-Fc)^Phh(k,mmin) + l/riE. Here, we use "in- 
trinsic" to mean the value of these functions when the effect 
of damping wing absorption from the HII regions is not in- 
cluded. Note that nh(mmin) is the number density of halos 
with m > mmin, i^hh(k, mmin) is the power spectrum of ha- 
los that have m > mmin (with the shot noise component 
subtracted off), and this expression for Pint,E assumes that 
there is no clustering of the foreground interlopers^ 

We ignore peculiar velocities in our analysis. Peculiar 
velocities most significantly affect Fourier modes oriented 
along the line of sight. In narrow-band LAE surveys, the 
majority of the signal is instead contributed by transverse. 



^ Even though this assumption might not be true in detail, the 
majority of observed emitters reside in halos of roughly the same 
mass in our model because current observations probe the tail of 
the luminosity function. 

^ For narrow band LAE surveys, foreground interlopers come 
from a limited number of redshifts because they must have a 
strong emission line that falls within the narrow frequency band. 
As a consequence of this redshift selection, the clustering of fore- 
ground interlopers is more significant than if the interlopers came 
from all intervening redshifts. However, unless the contamination 
fraction of a LAE survey is unusually high, the extra clustering 
arising from the foreground interlopers should be much smaller 
than the contribution from the LAEs, justifying the no-clustering 
approximation. 
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long wavelength modes. Moreover, in linear theory the pecu- 
Uar velocity terms that contribute to Pint,E are suppressed 
by a factor ^ 1/6, where b is the n^;- weighted linear bias 
of emitters. Based on the abundance of observed z = 6.6 
LAEs, a reasonable estimate is 6 7 in present surveys. 

This model gives the number density and clustering of 
emitters with m > mmin. However, surveys do not directly 
infer the mass of an emitter but instead measure its Lya 
luminosity. In the absence of dust, we can map the ionizing 
luminosity from a so urce in our simula tions to its luminosity 
in Lya photons via (|Osterbrocklll989l ) 

=0.67^^(1 -/esc) iVion, (3) 

where = 10.2 eV, /esc is the fraction of ionizing photons 
that escape into the IGM, and iVion is the production rate 
of ionizing photons. Assuming that equation (|3]) holds, that 
/esc = 0.02, and that Lobs,E = 0.5 Lq, (where the factor 
0.5 roughly accounts for resonant absorption in the IGM), 
the observed Lya luminosity of the sources in our fiducial 
reionization model [model (i) in the next section] is Lobs,E = 
5 X 10^^m/(10^^ Mq) erg s"^ The emitters that have been 
observed at z = 6.6 h ave Lobs,E = (2 — 20) x lO'^^ erg s~^ 
(|Ta niguchi et al. 2005). Lowering or raising /esc - which is 
currently unconstrained - can decrease or increase Lobs,E for 
a halo of fixed m. 



photon sinks. All minihalos with m > 10^ Mq absorb in- 
cident ionizing photons out to their virial radius until they 
are photo-evaporated. This absorption cross section is larger 
than the cross section found in ra diative- hydrodyn amic sim- 
ulations of minihalo evaporation (|lliev et al.ll2005l V suggest- 
ing that this model overestimates the impact of minihalo 
absorption. The photo-evaporation timescale is roughly the 
sound- crossing timescale of a halo0 

Figure [1] displays slices through simulations using mod- 
els (i), (ii), and (iii). The white regions are ionized and the 
black are neutral. Model (ii) results in the largest HII re- 
gions because it has the most biased sources, whereas model 
(iii) produces the smallest bubbles, with the maximum bub- 
ble radius restricted to be roughly the mean free path for 
ionizing photons to intersect a minihalo. 

One piece of physics that is missing from these three 
reionization models and that has has not been quantified in 
simulations of reionization is the effect of a duty cycle for the 
ionizing sources on the morphology. It is probable that the 
galaxies form massive stars and contribute ionizing photons 
only during certain periods. In Appendix[Cl we demonstrate 
that the duty cycle of the sources does not affect the mor- 
phology of the ionized regions for most realistic reionization 
models. 



2.4 Three Reionization Models 

We consider three models for the reionization history. These 
models bracket the possible range of morphologies for reion- 
izatio n by POPII-like stars (as shown in .McQuinn et al.l 
l2006l ): 

(i) All halos above mcooi contribute ionizing photons at a 
rate that is proportional to their mass m. The ionizing lumi- 
nosity of source halos is iVion(m) = 3 x 10^^ [m/(10^° Mq)] 
ionizing photons s~^. The scaling iVion m assumes that 
the star formation rate is proportional to the amount of 
gas within a galaxy. The normalization of A^ion is chosen 
such that reionization ends at z — 7. Given the uncertainty 
in /esc and the star formation rate in high-redshift galaxies, 
the range of possible normalizations is vast. Fortunately, the 
morphology of reionization dep ends very weakly on th e nor- 
malization of iVion , as shown in iMcQuinn et al.l (|2006l l . 

(ii) Halos more massive than mcooi contribute to the 
production of ionizing photons, with the ionizing luminos- 
ity of the sources scaling as halo mass to the 5/3 power 
(i.e., more massive halos dominate the production of ion- 
izing photons than in model (i)). This scaling is chosen to 
match the relationship between star formation efficiency and 
galaxy m ass that is observed in low-redshift dwarf galaxies 
(|Kauffmann et al.l 120031 ) as well as the star formation rate 
found in theoretical studies that include superno va feedback 
(jPekel Wooll2003l : ISD7ingel Hernauistll2003l ). The total 
budget of ionizing photons released is calibrated such that 
reionization ends at 2; ~ 7 as in model (i). 

(iii) Absorption by minihalos shapes the morphology of 
reionization. We use the same source parametrization as in 
model (i) except that the sources are twice as luminous so 
that reionization ends at z ~ 7. While minihalos do not 
contribute ionizing photons in this model, they do act as 



3 LAE MAPS 

Figure [2] shows mock LAE surveys created using snapshots 
from the simulation of model (i) in the 94 Mpc box and with 
a depth of 130 A or approximately 35 Mpc. The dimensions 
of these mock surve ys are roughly the same as the z — 6.6 
Suba ru Deep Field (^Taniguchi et al. I l2005l :l Kashikawa et al. 
"2006'). Each panel would subtend 0.6 degrees or roughly 
the solid angle of the moon. These LAE maps are gen- 
erated using Method 2 in Section 12.21 and assuming that 
25% of halos host emitter^ (i.e., /e = 0.25) and that 
only those Lya emitting halos can be observed that have 
m exp(— ra(z^o)) > 7 x 10^° Mq. This model corresponds to 
setting mmin = 7 x 10^° Mq. Note that Lint e cx: m in model 

The top panels in Figure [2] are the projected ionization 
maps for Xi = 0.3 (left panel), 0.5 {middle panel), and 0.7 
{right panel) at redshifts z = 8.2, 7.7, and 7.3, respectively. 
The projected map is completely ionized in white areas and 
neutral in black ones. The middle row of panels shows the 
intrinsic distribution of emitters (or what would be observed 
if Xi ^ 1), and the bottom row of panels shows the observed 
distribution. (Compare the bottom panels with their corre- 
sponding top panels to see the reionization-induced modu- 
lation.) The "Intrinsic" ue in Figure [2] is a few times higher 

^ We use the fitting formula in llliev et al] (|2005l ) for the evapo- 
ration timescale, which parameterizes the evaporation timescale 
as a function of redshift, halo mass, and incident ionizing flux. 
^ Twenty-five percent is chosen to match the fraction of Lyman- 
break galaxies at 2; = 3 that meet the se lection criteria to b e 
detected in a narrow band survey as a LAE (Shapl ev et aIll2QQ3h . 
although it is probable that Lyman break galaxies have some duty 
cycle as well and that the fraction of galaxies that emit in Lya 
evolves with redshift. 
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Model (i) Model (ii) Model (iii) 







Figure 1. Slices through the middle of the 186 Mpc simulation box for the three models. The top row is at Xi = 0.3, the middle is 
at Xi = 0.5, and the bottom is at Xi = 0.8. Model (iii), in which minihalos limit the photon mean free path, has the smallest bubbles, 
whereas model (ii), in which the sources are the most biased, has the largest. 



then the tie of the z = 6.6 Subaru De ep Field (SDF) photo- 
metric sample (|Taniguchi et al. The "Observed" tie 
when Xi = 0.5 in Figure [2] is comparable to the ue of the 
z = 6.6 SDF photometric sample. Note the higher degree of 
clustering in the "Observed" panels relative to the "Intrin- 
sic" panels. We investigate the detect ability of this clustering 
in Section [5] 

Figure [3] is the same as Figure [2l but for a futur- 
istic survey that is sensitive to Lya emitting halos with 
m exp(— Ta(z^o)) > 1 X 10^° Mq. There are approximately 
2000 sources in the "Intrinsic" panels, but there can be sig- 
nificantly fewer in the "Observed" panels. The large-scale 
modulation of emitters by the bubbles in Figure [3] allows 
observations to image the neutral holes in the map with 
the LAEs. If the Lya luminosity scales roughly with m, as 
it does for the sources in these calculations, the mock sur- 



vey depicted in Figure [3] would entail a ^ 50 times longer 
observation than that of Figure [2] (assuming that the obser- 
vation is photon- limited). A space mission with a wide field 
of view (FoV) camera could provide such a deep survey in 
a much shorter observation. The James Webb Space Tele- 
scope (JWST), set to be launched in 2013, will be able to 
easily provide such a deep image. However, the current high- 
redshift program for JWST specifies a deep survey that will 
cover less than one-tent h of the solid angle s ubtended by a 
single panel in Figure [3l ([Gardner et al.|[2006l V With such a 
small field, it will be difficult for JWST to probe the HII 
regions during reionization. The small square in the lower 
left-hand panel represents the FoV of JWST. 
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Figure 2. Top panels show the projection of Xi in the survey volume. In the white regions the projection is fully ionized and in black 
it is neutral. The left, middle, and right panels are for z = 8.2 (xi = 0.3), z = 7.7 (xi = 0.5), and z = 7.3 (xi = 0.7). The middle and 
bottom rows are the intrinsic and observed LAE maps, respectively, for /e = 0.25 and assuming that we can observe unobscured emitters 
with m exp(-rc(i/o)) > 7 x lO^^M©. (Note that 

-^int.E m.) The observed distribution of emitters is modulated by the location of the 
HII regions (compare bottom panels with corresponding top panels ). Each panel is 9 4Mpc across (or 0.6 degrees on the sky), roughly 
the area of the current Subaru Deep Field (SDF) at z = 6.6 (Kas hikawa et all 120061 ). The depth of each panel is AA = 130 A, which 
matches the FWHM of the Subaru 9210 A narrow band filter. The number densities of LAEs for the panels in the middle row are few 
times larger than the number density in the SDF photometric sample of 2; = 6.6 LAEs. 



4 THE LUMINOSITY FUNCTION 

The luminosity function of LAEs depends sensitively on the 
morphology of HII regions during reionization. If all the HII 
regions are smaller than 1 pMpc (such that Ta(z^o) > 1)), 
only a small fraction of LAEs would be observed compared 
to the number that would be observed if the Universe were 
fully ionized. In fact, even if the bubbles are a few times 
larger than 1 pMpc, many emitters will be obscured. This is 
because the steep, decreasing nature of the luminosity func- 



tion means that the majority of emitters have luminosities 
that are within a factor 2 of the detection threshold, requir- 
ing for many LAEs that ra(z^o) be significantly less than 
unity to be observed. 

As reionization proceeds, larger HII regions will form, 
allowing more LAEs to appear out of the dark. An extremely 
rapid decrease in tie would be difficult to attribute to evo- 
lution in the intrinsic properties of the LAEs rather than 
to reionization. Since the LAE luminosity function has been 
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Figure 3. Same as Figure [21 but for a futuristic LAE survey that can detect halos down to m exp(— Tq, (z^o)) > 1 x 10-*^^ (note 
Lint,^ cx: m) and assuming /e = 0.25. The large-scale modulation of LAE by the HII bubbles is clearly apparent in this survey. The 
square in the lower left-hand panel represents the 3^ X 3^ FOV of JWST drawn to scale. 



measured at z = 6.6 (and tentatively at z = 7) and v^ill be 
constrained at even higher redshifts in the coming years, it 
is important to understand the signature of reionization in 
the luminosity function. 



The current data on LAE luminosity functions may 
indicate that re i oniza tion is happening at z — 6.6. 
iKashikawa et al.l (120061 ) finds that there is a suppression 
in the bright end of their measured luminosity function at 
z — 6.6 relative to that dX z — 5.7 at 2-cr significance. 
See the thick solid curve in Figure [5] for the ratio of the 
z — 5.7 and z — 6.6 best-fit LAE luminosity functions 
along v^ith an estimate for the 1-cr shot noise errors on 



this ratio (IShimasaku et aP l2006l : iKashikav^^a et al.l ^ 
IKashikawa et al.l (|2006') propose that the suppression of the 
high luminosity end may imp ly a chang j e in t he ionization 
state of the IGM. In addition, live et al I ([2007l ) finds an ad- 
ditional factor of few suppression in the luminosity function 
at ^ = 7 for emiter s wit h Lobs,E > 1 x lO'^^ erg s~^. 

Dijkstra et al. (2006) suggest that there is a more mun- 
dane explanation for this suppression - the evolution of the 
halo mass function. Employing a similar model for LAEs 



^ The cosmic variance errors are highly correlated between dif- 
ferent luminosity bins, and a 1-cr cosmic fluctuation will raise or 
suppress this ratio by ~ 50% for an ionized universe (see Section 



[6]). 
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Figure 4. Number density of LAEs with Lo^g e > L dit z = 
6.6 and for several different volume- averaged ionization fractions. 
These curves are calculated from a simulation of model (i) in 
the 186 Mpc box. The top axis shows the number density with 
m exp(— rcK(z/o)) > ^min- The mapping between m and Lobs,E 
is discussed in the text. 



to ours, they argue that the cosmo logical evolution of the 
halo mass function between z = 5.7 and z = 6.6 can ac- 
count for the observed suppression a nd that it is unneces - 
sary to invoke reionization. However, if lPiikstra et al.l (|2006l ) 
were to include dispersion in the luminosity for a given halo 
mass, the high mass end of the luminosity function would 
be less dependent on the evolution of the halo mass func- 
tion, possibly altering this conclusion. Also, the model in 
iDiikstra et all (|2006l ) favors /e - 1 in order to fit the data 
at z — 6.6 whereas they find that /e ^ 1 provides a bet- 
ter fit at z = 5 .7. A similar trend for the evolution of /e 
is also found in IStark et"al ] £2007), even though such dra- 
matic evolution in the intrinsic properties of emitters is not 
expected. 

For most of the calculations in this section, we as- 
sume the simple model for LAEs discussed in Section 
12.31 Halos that are active LAEs have ivint,E(^) = 5 x 
lO^^(/esc/O.O2)-^[m/(lO^^M0)] erg s"^ in model (i). We 
set /esc = 0.02 to roughly match the observed abundance of 
LAEs for /e = 0.1 and mmin = 5 x 10^° M©. For the other 
two models, we take Lint,E(^) to be the same function as 
in model (i). To be consistent with equation JS]) and the 
A^ion used in models (ii) and (iii), fixing Lint,E(^) requires 
adjusting /esc slightly. 

The curves in Figure |4] represent the observed number 
density of LAEs with luminosities that satisfy Lobs,E > L. 
These curves are calculated from the 186 Mpc simulation 
of model (i) at z = 6.60 During reionization, the luminos- 
ity function is suppressed by a factor that does not depend 
strongly on L. Figure [5] plots the ratio of the luminosity 
function at various Xi for the three models. As with model 



^ To perform this calculation, we use the halo field at z = 6.6 in 
the 186 Mpc box, but we use the ionization field from the simula- 
tions of model (i), which are generally at slightly higher redshifts. 
This is justified because the ionization field, when comparing at 
fixed Xi, is essentially independent of the redshift where this Xi 
is reached (McQuinn et al. 2006. V 
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Figure 5. Ratio of the number density of LAEs with I/obs,E > L 
for different x^, nE(> L,Xi), to the number density for Xi ~ 1, 
nE(> L, 1). The thick solid curve with the 1-a errors in the top 
panel is the ratio of the z = 6.6 and 2; = 5.7 SDF luminosity func- 
tions. The effect of reionization in all models is approximately a 
uniform suppression of the luminosity function that is indepen- 
dent of luminosity. Model (ii) has the largest bubbles such that 
the luminosity function is the least suppressed at fixed Xi, whereas 
model (iii) has the smallest bubbles such that the luminosity func- 
tion is the most suppressed. 
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Figure 6. The same as Figure [S] but the curves are computed 
assuming model (i). Top Panel: Curves compare the suppression 
in the fiducial model in which L (x m (source model A) to a model 
in which L is independent of m (source model B). Bottom Panel: 
Curves labelled "wind" assume = 400 km s"-*^ and labeled 
"infall" assume -Uw = '^^vir- 



(i), the luminosity function in the other models is suppressed 
by a largely luminosity-independent factor at each Xi. The 
sharp cutoff at high luminosities in this ratio owes to the 
finite simulation volume. Model (ii) has the largest bubbles 
so that the luminosity function is the least suppressed at 
fixed Xi, whereas model (iii) has the smallest bubbles such 
that the luminosity function is the most suppressed. 

We also calculated the suppression of the luminosity 
function for a wind model with — 400 km s~^. This value 
for is motived by the measured average ve locity offset of 
the Lya line in strong emitters at ^ ~ 3 (Shap lev et al.l 
I2OO3I I. This model probably yields an upper limit for the 
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Figure 7. Probability of having Ta{h'o) < r. The thin curves are 
for LAEs with 1 x 10^° Mq < m < 2 x 10^° M©, and the thick 
curves are for LAEs with 5 x 10^° Mq < m < 1 x 10^^ Mq. 
The sohd curves are for Xi = 0.3, the dotted curves are for Xi = 
0.5, and the dot-dashed curves are for Xi = 0.7. Assuming the 
observational threshold of m > 1 X 10^^ exp(— (z/q)) only 
the emitters that contribute to the portion of the thin (thick) 
curves left-ward of the thin (thick) vertical lines are observed. 



effect of winds at z > 6, where the galaxies are much less 
massive than at z ^ 3 and, therefore, not able to power 
such strong winds. A wind causes a redshift of the Lya Une 
because Lya photons lose energy when they scatter off the 
baryons in the wind. Winds make the bubble size needed for 
a fixed damping wing optical depth shrink. See the curves 
labeled "wind" in the bottom panel in Figure [6]and compare 
with the other curves with the same Xi. This comparison il- 
lustrates the effect of a winds on the luminosity function 
can be significant for x^ < 0.3, but tend to be unimpor- 
tant when Xi > 0.6. The effect of winds on the luminosity 
function decreases with increasing Xi because as the bub- 
bles grow (increasing Xi) the ratio of the wind velocity to 
the Hubble flow velocity at the bubble edge, v^^/ {H{zg) Rt), 
decreases. 

We also consider a model in Figure [6] where we set 
'^^w = Vvir{m), where Vvir{m) is the halo circular velocity at 
the virial radius. This model is meant to emulate analytic 
calculations in which resonant absorptions of infalling mate- 
rial can obscure all wavelengths bl ueward of A ( 1 + ^;vir (m) / c) 
(|Santosll200j : IPiikstra et al]l2007l l . As with winds, infall has 
a relatively minor effect at Xi > 0.6, but can reduce the sup- 
pression of the luminosity function for Xi < 0.3 (Figure [6]). 

The uniform suppression of the luminosity funct i on wa s 
also found in the analytic st udies of iFurlanetto et aD (|2006l ) . 
The explanation provided in Furlanetto et aL ( 20061 ) for this 
uniform suppression is that the most massive (most biased) 
sources sit in the largest bubbles, which results in the least 
attenuation of their Lya line. In addition, the most mas- 
sive LAEs are the most luminous (at least in our model), 
requiring a larger than an average mass halo does to be 
obscured. However, at the bright end of the luminosity func- 
tion, even a slight decrease in the average luminosity of the 
LAEs causes the luminosity function to change rapidly be- 
cause of its steep slope. These effects sum to give roughly 
the same suppression at the bright end as at the faint end. 
In what follows, we investigate this explanation in detail and 



quantify the effect of our assumptions on our predictions for 
the evolution of the luminosity function. 

Figure [7] illustrates the dependence of damping wing 
absorption on halo mass, plotting the cumulative probability 
distribution of Tc(z^o), P{Ta{i^o) < r). The thick curves are 
P(Ta(z/o) < t) for LAEs with 5 x 10^° M© < m < 1 x 
10^^ M©, and the thin ones are P(Ta(z/o) < r) for LAEs with 
1 X 10^° Mq <m <2x 10^° Mq. This figure demonstrates 
that the most massive LAEs have smaller because they 
sit in larger bubbles than less massive LAEs. This effect 
decreases with increasing Xi because the bubbles become 
less biased as they grow. (Compare the difference between 
the two soUd curves - = 0.3 - to the difference between 
the two dot-dashed curves - = 0.7.) In addition, the most 
massive LAEs are the most luminous in our model, which 
requires a larger than less massive LAEs do for their Lya 
luminosity to fall below the survey threshold Lmin- 

Figure [71 suggests that our predictions for the luminos- 
ity function depend on the scaling of the LAE bias with 
luminosity. This is worrisome because this scaling is very 
uncertain. How dependent are our predictions on the sim- 
ple mapping we assume between halo mass and Lya lumi- 
nosity? One plausible extreme LAE model is, rather than 
i^int,E oc m, for Lint,E to be independent of halo mass. To 
achieve this, we kept the same halo positions but random- 
ized the luminosities of the halos while maintaining the same 
luminosity function as in the fiducial Lya emitter model. We 
refer to this source model as model B and the fiducial model 
as model A. The top panel in Figure [6] compares the sup- 
pression of the luminosity function in these two models. The 
suppression in the two models is similar, particularly for the 
Xi = 0.62 case. The agreement between these two models 
suggest that our predictions are not strongly dependent on 
our prescription for the luminosity function. 

When Xi — 0.33, the luminosity function is more sup- 
pressed at the high mass end for source model B than for 
source model A, imparting some scale dependence in the 
suppression for model B. The scale dependence induced by 
model B can be thought of as the maximum scale depen- 
dence that can be imparted by reionization since bias is 
uncorrelated with luminosity in this model. Therefore, it 
is difficult for reionization to be solely responsible for the 
scale-dependent suppression that may have been observed 
in Kashikawa et al. (2006) E3 

If a LAE survey suffers from significant contamination 
or incompleteness, these systematics will affect the normal- 
ization and shape of the luminosity function, complicating 
any inference that the evolution is due to neutral regions. 
Also, intrinsic evolution in the source properties with red- 
shift could be difficult to distinguish from reionization. If 
evolution in the observed luminosity function indeed owes 
to reionization, then reionization will also increase the clus- 
tering of observed emitters, whereas systematic effects and 
intrinsic evolution cannot change the clustering in the same 



It is possible to induce more of a scale dependence in the sup- 
pression by making the luminosity function steeper than we have 
assumed at large luminosities. However, we find that in practice 
a steeper bright-end luminosity function does not result in much 
additional suppression. This result owes to the large dispersion in 
T(x for halos of fixed Lobs,E- 



© 2006 RAS, MNRAS 000, 1 



manner. The impact of reionization on clustering is dis- 
cussed in the following section. 



5 CLUSTERING 

Suppose that upcoming observations confirm the substan- 
tial decrease in bet ween z = 5.7 and z = 6.6 found in 
Kashikawa et al. 1 (|2006h . This evolution could be explained 
in three ways (or some combination thereof): 

(i) Decreasing Duty Cycle: The number density of halos 
that host emitters is decreasing, but the observed emitters 
sit in halos of the same mass at z = 6.6 as at z = 5.7. 

(ii) Increasing Halo Mass: The average halo mass of LAEs 
is increasing between z = 5.7 and z = 6.6. 

(iii) Reionization: Patchy reionization is still occurring at 
z = 6.6. Neutral regions are obstructing the line of sight to 
some emitters, decreasing their observed abundance. 

In case (i), the clustering properties of the sources will be 
essentially unchanged in the 150 million years between z — 
5.7 and z = 6.6. Only in cases (ii) and (iii) can the amount 
of clustering increase significantly. If a high-redshift LAE 
survey is able to distinguish case (ii) from case (iii), it will 
be capable of determining whether reionization is happening 
at z = 6.6. 

The curves in Figure [8] are the 3-D power spectrum of 
Se = riE/fiE — I, calculated from a simulation of model (i). 
The shot noise component of this power spectrum has been 
removed from these curves. Note that these curves depend 
only on mmin and the ionization field; they do not depend 
on /e. We plot the 3-D rather than the 2-D power spectrum 
because the 3-D power spectrum makes use of all the two- 
point information that is available in our simulation volume, 
minimizing cosmic variance. Note that all of our conclusions 
about the effect of reionization on the LAEs would be the 
same if we considered the 2-D power spectrum. The curves 
in the middle panel in Figure [8] represent the case that the 
observed LAEs have m exp(-rc(z^o)) > 7 x 10^° M©. (Note 
that Lint,E ^ m.) This threshold along with /e ~ 0.1 yields 
the riE measured in the SDF z — 6.6. 

The curves in the bottom panel in Figure [8]represent the 
case in which LAEs are observed in less massive halos than 
the middle panel, halos with m exp(— rc(z^o)) > Ix 10^° M©. 
The curves correspond to the signal for surveys that are 
more sensitive than Subaru such as JWST (unless /e < 
0.01). The LAE power spectrum at fixed Xi has a slightly 
lower amplitude when computed from more abundant halos. 

During reionization, the power spectrum of LAE fluc- 
tuations changes rapidly. Between Xi ^ 1 (when z = 6.7 in 
the simulation of model (i)) and Xi — 0.5 (when z = 7.5), 
the amplitude of the power spectrum increases by roughly a 
factor of 3 — 4 (see middle panel in Fig. (8]). The amplitude 
increases by another factor of ^ 2 by = 0.3 (when z = 8). 
If LAE surveys detect a rapid increase in the amplitude of 
the power spectrum with redshift, it would be difficult to 
attribute this to anything other than reionization. 

The straight, solid diagonal lines in Figure [8] represent 
the shot-noise power spectrum for tie = 1 x 10~^ Mpc~^ 
and nE = 5 X 10~^ Mpc"^. The former riE IS approxi- 
mately the riE measured in the SDF photometric sample 
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Figure 8. Dimensionless 3-D power spectrum of = n^/hE — l, 
calculated using the 186 Mpc simulation of model (i). Bottom 
Panel: The power spectrum of emitters with m exp(— Tq, (z^o)) > 
1 X IQ-*^^ ^0- All curves are from redshifts between z = 6.9 and 
z = 8.3. Middle Panel: The power spectrum of emitters with 
m exp(-ra(z^o)) > 7 x lO^o Mq, and fiint,E ~ 10~^ /e Mpc~2. 
Top Panel: The intrinsic clustering of emitters (or, equivalent ly, 
the clustering when Xi ^ 1). The thick curves are if emitters can 
be observed in halos with m > 1 x IQ-'^-'^ Mq at z = 6.9 (dashed 
curve, fiE = 1 X 10~^ /e Mpc"*^) and at z = 8.0 {solid curve, 
fiE = 5 X 10~^ /e Mpc~^). The thin solid and dashed curves 
are the same but in halos with m > 1 x 10-*^^ Mq (with tt-e = 
0.03/eMpc-3 at z = 6.9 and nint,E = 0.02 /e Mpc-^ at z = 8, 
respectively). All Panels: The straight, vertical line corresponds 
to k = 27T h , where R = 80 Mpc - roughly the angular extent 
of the z = 6.6 LAE survey in the SDF. The diagonal solid lines 
are the shot-noise power spectrum for tt-e = 1 x 10~^ Mpc""^ and 
fiE = 5 X 10~^ Mpc~^. The former is approximately the value 
of fiE in the SDF photometric sample. 

(jKashikawa et al.ll2006h . Notice that for a survey with ei- 
ther value of nE, clustering can be detected on large scales 
provided that the survey volume is large enough (i.e., the 
shot-noise line is below the other curves). The vertical 
dashed lines in Figure [8] correspond to /c = 2ti hR~^, where 
= 80 Mpc - roughly the angular extent of the z — 6.6 
LAE survey in the SDF. These curves imply that if reion- 
ization is happening dX z — 6.6, fluctuations in the LAE field 
can be imaged even in current programs. 

The results in the bottom two panels of Figure [8] should 
be contrasted with what observations would see if, instead 
of the increased clustering owing to reionization, it was en- 
hanced by evolution in the intrinsic properties of the LAEs. 
A similar evolution in the clustering properties to the evolu- 
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Figure 9. Top Panel: Dimensionless 3-D power spectrum of 
(5e for three different models of reionization, calculated from 
snapshots that have Xi = 0.5 and assuming m exp(— rQ;(z/o)) > 
1 X 10"'^'-' Mq. In addition, a curve for the intrinsic (or Xi ~ 1) 
power spectrum is included for comparison. Bottom Panel: Com- 
parison of the Xi ~ 0.5 LAE power spectrum of source model A 
(in which Lint,E ~ ^) to both the LAE power spectrum of a wind 
model that uses source model A and = 400 km s"-*^ and to the 
LAE power spectrum of source model B (in which I/int,E is inde- 
pendent of m). The curve labeled "intrinsic" uses source model 
A, and all curves in the bottom panel use reionization model (i). 

tion caused by reionization could conceivably be produced 
if the bias of the sources were to increase with redshift - i.e. 
if the LAEs sit in more massive halos at z = 6.6 than at 
z = 5.7. 

In the top panel in Figure [H the thin curves have a 10 
times smaller mmin than do the thick ones. For fixed /e, 
there are over 100 times more emitters in the survey from 
which the thin curves are calculated compared to the survey 
from which the thick curves are computed. A change m TTlmin 
that leads to a change in tie by a factor of over 100 yields a 
smaller variation in Pe than the change in Pe between when 
Xi ^ 1 and when Xi = 0.5, even though tie differs by only 
a factor of 3 between these curves (fixing /e). Reionization 
causes the clustering of LAEs to evolve much more quickly 
than is possible with intrinsic evolution. 

Let us develop a general understanding of the magni- 
tude by which intrinsic evolution of the LAEs can change 
the amount of clustering. On large scales, Pe = ^^Pdm, 
where b is the intrinsic large-scale bias and Pdm is the 
dark matter power spectrum. The large-scale bias calcu- 
lated using Press- Schechter theory at z = 7 is 6^ = 28 for 
mmin = 1 X 10^° Mo, = 58 for mmin = 5 X 10^° M©, 
b^ = 65 for mmin = 7 x 10^° M©, and b^ = 77 for 
^min = 1 X lO"*^"*^ Mq. The ratio of the Pe for ~ 1 and for 
different mmin is typically smaller than the ratio between Pe 
for models in which Xi differs by Axi ^ 0.3 . 

The top panel in Figure [9] compares the power spectrum 
of for reionization models (i), (ii), and (iii) and for Xi = 
0.5 with mmin = 1 x 10^° M©. Model (iii) has the least power 
on large scales - despite having the fewest observed emitters 
of the three models - because it has the smallest bubbles. 
However, the differences between the curves for = 0.5 are 
not as substantial as the differences between Xi = 0.3, Xi = 
0.5, and Xi ~ 1. Therefore, it will be easier to constrain Xi 



to ^ 30% with emitter clustering than it will be to constrain 
the details of the reionization process. 

The bottom panel in Figure [9] compares the power spec- 
trum of & computed from a LAE field that uses the fiducial 
source model (model A) and from a LAE field that uses 
source model B, which was discussed in Section |4l where 
luminosity is independent of halo mass for all halos above 
^min- All curves are calculated with mmin = 1 x 10^° Mq. 
The bottom panel in Figure [9] also plots the power spectrum 
for a LAE model with — 400 km s~^ for all emitters. 
Winds suppress the amount of clustering. 

If a LAE survey detects excess power on large scales, 
a skeptic might contend that this is caused by a rare large- 
scale structure in the survey field rather than by reioniza- 
tion. Fortunately, there is a simple test that may help to 
distinguish between these two hypotheses: Survey the field 
using a different selection criterion, such as Ha emission or 
by the Lyman-break technique. If the reionization hypothe- 
sis is correct, less clustering is expected in the second survey, 
and if we cross correlate with the Ly<^ survey, again the same 
excess clustering will be present 

This test does not rule out the possibility that the 
LAEs are intrinsically more biased than the galaxies se- 
lected with the other selection criterion. This possibility 
may be surprising because at lower redshifts LAEs are as- 
sociated with younger, less biased systems than the Lyman 
break galaxies. However, let us suppose that this circum- 
stance is the case. On scales where the intrinsic bias is lin- 
ear and shot noise is unimportant, the Fourier space fiuctu- 
ations of the two galaxy fields should have the same phase 
at each k if ^ 1. If this is not true, then this would 
be evidence for reionization. Of course, on large enough 
scales the emitter field, &, will have the same phase as the 
galaxy field, even during reionization, but the presence 
of HII bubbles will extend the range of scales over which 
the phases do not agree. The cross correlation coefficient 
r{k) = compares the phases 

of two fields. If Xi — 0.5 and mmin = 5 x 10^° M©, we find in 
model (i) that r ~ 0.9 at A: = 0.1 Mpc~^, whereas if ~ 1 
then r^latA: = 0.1 Mpc~^ (even if the mmin differ by a 
factor of a few between the and 5g fields). 

Another potential test to check whether the fiuctua- 
tions in the LAE field owe in part to reionization is to note 
that the HII regions break the rotational isotropy of the LAE 
field because Lya absorption depends on the ionization state 
of the IGM in only the line-of-sight direction We have 



There is a second-order effect owing to Jeans mass suppres- 
sion in ionized regions inhibiting galaxies with m < 10^ Mq from 
forming. A high redshift galaxy would be composed of many of 
these smaller mass galaxies owing to past merger events, and 
therefore its Lyo; emission mig ht be influenced by the local reion- 
ization history (Babich &: Loebll2006l : iPritchard et al.ll2007l ). This 
effect may result in reionization also influencing the clustering of 
all galaxies, not just LAEs. Quantifying the relevance of this effect 
is difficult, but it is certainly much smaller than the modulation 
we consider. 

^'^ Physically, a uniform distribution of LAEs and a spherical HII 
region will result in a non-spherical observed distribution of LAEs 
within the HII region, with the distribution of LAEs truncated at 
~ 1 pMpc from the edge of the HII region along the line of sight 
direction. In reality, the distribution of emitters and the shape of 
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Figure 10. This plot investigates the effects that shape 
kP' P^{k) / 2tt'^ . The dashed curve is the intrinsic LAE power 
spectrum for mmin = 1 x 10-'^'-' Mq. The thick sohd curve is 
k^P^{k) /2tt'^ for Xi = 0.5 and mmin = 1 X 10^° Mq. The dotted 
curve is the same as the thick sohd except that each of the halos 
has been randomly displaced in the box prior to computing the 
observed LAE field. The dot-dashed curve is the same as the thick 
solid except that the halo field is uniformly translated by some 
arbitrary distance relative to the ionization field (to make Xi and 
6g uncorrelated) prior to computing the observed LAE field. The 
shot noise contribution has been subtracted from the curves that 
represent Pe- The thin solid curve is five times the power spec- 
trum of Sx^ , included to compare its spectral dependence with 
the other curves. 

investigated the statistic (|fe(k • n)|^) for different orienta- 
tions of fi and have found no significant signature of angular 
anisotropy in our simulation cube (neglecting redshift-space 
distortions owing to the peculiar velocity field). 

In this section, we focused on the riE-weighted 3-D 
power spectrum. The reader might wonder why we did not 
use the luminosity-weighted power spectrum (or even some 
more general weighting) instead, which can contain more 
information about reionization. We find that a luminosity- 
weighted power spectrum results in a similar spectrum of 
fluctuations. The reason for this is that most emitters in a 
survey are near the detection threshold because of the de- 
creasing nature of the mass function. It only takes minor 
attenuation for most LAEs to not be observed, such that 
the major source of fluctuations derives from whether an 
emitter is detected or not. 

5.1 Factors that Shape the LAE Power Spectrum 

Let us develop a simple model for the LAEs to understand 
the reason that reionization has such a large effect on Pe- 
We will not use this model to predict Pe, but rather to 
provide a framework with which to understand the different 
effects that contribute to Pe during reionization. For this 
model, we assume that only emitters in bubbles of radius 
/ are observed, where / ^ 1 pMpc. In reality, / depends on 
the emitter luminosity and on the minimum luminosity of 
the survey. We can write the observed number density of 

HIT regions is more complicated than in this spherical example, 
which reduces the magnitude of this effect. 
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LAEs at position x as nE(x) = rig{x.) x7(x), where Xi is the 
ionization field smoothed over a sphere of radius /, and Ug 
is the local number density in galaxies that emit in Lya. If 
a bubble is defined as a fully ionized sphere then 7 = oo. 
However, in practice regions are never fully ionized, and for 
our purposes it suffices to leave 7 as a free parameter. 

In this model, the correlation function of LAEs can be 
written as (omitting constants) 



(5e 5e'> = X' [{x] x]') + 2 {x]S, x]') + {x]S, x]'d/}) , 

(4) 

where & and 6g are the overdensities in observed LAEs and 
Lya emitting galaxies, respectively, and where X = fig/n^. 
Even though this model is simplistic, if we had the full model 
for this correlation function it would have a similar decom- 
position to the decomposition seen in the right-hand side 
(RHS) of equation (g)). 

Which terms on the RHS of equation (|4]) shape &? The 
first term is generated only by the bubbles. On scales where 
this term is dominant, Pe is independent of the intrinsic 
clustering of the LAEs. To investigate the importance of this 
term, we randomly displaced all the halos in the simulation 
box such that the intrinsic halo field is Poissonian prior to 
computing the LAE field, fc(x). This operation makes the 
second and third terms in equation (UJ) zero for finite sepa- 
rations. The power spectrum of & using the displaced halo 
field for Xi = 0.5 is given by the dotted curve in Figure [10] 
(with the shot noise contribution subtracted off). Compare 
this curve to the thick solid curve, which represents Pe com- 
puted from the true halo field. On large scales, these curves 
agree fairly well, implying that the first term on the RHS of 
equation (|4| is important at tens of Mpc scales and greater. 
This result qua litatively agrees with the analytic study of 
iFurlanetto et al . (2006), which found that the bubbles dom- 
inate the large-scale emitter clustering. 

The thin solid curve in Figure [TOl is five times the power 
spectrum of dx^ (the normalization is chosen to facilitate a 
comparison of the shape of this curve with the shape of 
other curves). As expected, the spectral shape of the thin 
solid and dotted curves are similar. However, unlike the thin 
solid curve, the dotted curve is consistent with zero at /c > 
2 h Mpc~^. This difference is because sub-pMpc features in 
the bubbles do not affect the LAE field. 

We have quantified the importance of the first term, 
but what about the second and third terms on the RHS of 
equation (|4])? These terms depend in part on the covariance 
between 5g and Xi. To investigate the importance this co- 
variance, we translated the halo field in the simulation box 
relative to the ionization field such that 5g and Xi become 
uncorrelated. Pe increases at all k hy ^ 50% from this op- 
eration (compare dot-dashed curve with thick solid curve in 
Fig. [10]), and the shape of Pe is maintained. 

Since Pe is not drastically changed by correlations be- 
tween 5g and x^, this result motivates the assumption that 
5n and x% are uncorrelated to understand Pe- The second 
term on the RHS of equation @ is zero with this assump- 
tion and the third term becomes X'^ {6g 5g') {x] x] ). With 
these simplifications, the Fourier transform of equation (|4]) 
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becomes 



+ j d^k^P,,(k-kOP5757(kO , (5) 



where d^k = {2n)~^ dP\<i. 

At /c > 1 /iMpc-\ Pe is a constant times the in- 
trinsic power spectrum (see the curves in Fig. [8]). This 
Umit is easy to understand in this model. On these scales, 
P57 57 ^ and Pe ^ / d3k^P^^(k - k^ P.^^t (k^ ~ 
X'^ Pgg{k) j d^\<i' P6ci5:^{\<i'). The x7 field is primarily com- 
posed of zeros and ones, such that / d^k^ P5757 (k^) — 
(x^x^')(0) ^ {x^) = where (x^x^')(0) is the cor- 

relation function evaluated at zero separation. Therefore, 
equation (|5|) becomes Pe ~ X Pgg. At Xi = 0.5, uncorre- 
lated 5g and xi results in X 10 in our calculations, and 10 
yields the small-scale increase between Pgg and Pe in Figure 
flOl (compare dashed and dot-dashed curves) . 

If Pe oc Pgg on small scales then the proportionality 
constant must be X to give the correct shot noise term. In 
model (i), when Xi = 0.5 then X ^ 3, and when Xi = 0.3 
then X 8. Interestingly, these factors yield roughly the 
small scale power increase we see in Figure (8] Therefore, 
the relation Pe ~ X Pgg at /c > 1 Mpc~^ seems to hold in 
general and provides a consistency check for the reionization 
hypothesis, where X can potentially be estimated using a 
slightly lower redshift LAE survey to derive fig . 



5.2 Detectability 

iKashikawa et al. I (2006') measured the angular correlation 
function w^{r) from the SDF photometric sample of 58 emit- 
ters. Interestingly, they find that w^{r) is consistent with 
no clustering (see the connected circles in Fig. [12] for this 
correlation function). This result may allow us to put con- 
straints on at z = 6.6. In addition, the z = 6.6 LAE 
sample will increase in size by a factor of five in the coming 
year with the Sub aru/XMM-Newton Deep Survey (SXDS) 
(|Ouchi et al.ll2005l ). 

To proceed, we estimate the mean value of w^{r) as 
well as the covariance in w^{r) between radial bins, Ce, 
from many mock surveys. Given the dimensions of the sur- 
vey, we generate these mock catalogs in as many spatially 
independent volumes as our 186 Mpc simulation box allows, 
and we do this computation for rc(z^o) calculated along the 
i, J, and k directions. (For the SDF, 20 spatially indepen- 
dent surveys can fit into the simulation box, leading to 60 
mock catalogs. In practice, we construct many times more 
overlapping mock catalogs than this number to obtain all 
the information that is available from the simulation box.) 
In addition, we compute w^{r) and Ce from all the ha- 
los above mmin in each mock survey region, and then we 
subtract out the shot noise contribution to determine the 
cosmological part of the covariance matrix. This procedure 
takes advantage of the fact that the cosmological contribu- 
tion to w^{r) and C e does not depend on /e, allowing us to 
reduce the uncertainty in our estimates for these quantities. 

Figure [TTI plots predictions for the correlation function 
of LAEs at ^ = 6.6 and mmi n 7 x 10^ ° M©. The pes- 
simistic estimate in .Kashikawa et al.l (|2006h for the contam- 
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Figure 11. Angular correlation function of emitters at z = 
6.6, assuming that observed emitters reside in halos with 
m exp(— ra(z/o)) > 7 x 10 Mq. The curves in the top panel 
are calculated in the same volume and with the same number of 
emitters, 58, as the SDF photometric sample. The bottom two 
panels are in a volume a slightly larger volume than the upcom- 
ing 1 sq. deg. Subaru/XMM-Newton Deep Survey (SXDS), with 
250 emitters in the middle panel and with 190 in the bottom 
one. The thick error bars owe to shot noise, and the thin owe 
to shot noise plus cosmic variance. To calculate these errors, we 
conservatively assume Fc = 0.25 in the top two panels (Fc = in 
the bottom panel). Current surveys can potentially distinguish an 
ionized universe (the curves labeled "intrinsic") from a universe 
with Xi < 0.5. 



inatio n fraction of the SDF p hotometric sample is Fc = 0.27, 
and Kashikawa et ajj (|2006[ ) estimates that the contamina- 
tion is probably closer to Fc = 0.16. We set Pc = for the 
spectroscopic survey curves {bottom panel) and Fc = 0.25 
for the photometric surveys (middle and top panels). We 
also lower the number of emitters in the mock spectroscopic 
surveys by 1 — Pc, where Fc = 0.25 is the contamination 
fraction in the mock photometric sample. Foreground con- 
tamination will bias the estimate for the measurement of 
the correlation function by the factor (1 — Pc)^. Rather than 
plot biased curves for we in Figure [TT] we instead divide 
the Poisson errors by the appropriate factor to account for 
contamination. 

The top panel in Figure [TT] displays the average cor- 
relation function for several clustering models, generated 
in mock surveys with the same dimensions as the SDF 
(34^ X 27^ X 130 A) and with 58 emitters - the number of 
LAEs in the SDF photometric sample. The thick error bars 
in Figure [TT] account for shot noise and the thin error bars 
include both shot noise and cosmic variance. Note that the 
cosmic variance errors in the top panel are important, par- 
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Figure 12. Angular correlation function of emitters at z = 6.6 
measured from independent volumes in the 186 Mpc simulation 
box of model (i) for x-i ~ 1 {solid curves) and for Xi = 0.5 
{dashed curves). The curves in the top two panels are calcu- 
lated in the same volume and with the same number of objects, 
58, as the SDF photometric sample, and the curves in the bot- 
tom panel are calculated from 250 emitters in an area that is 
slightly larger than the upcoming SXDS. These curves assume 
that m exp(— rQi(z^o)) > 5 x 10^^ Mq, and there are 8 curves in 
each set. The thick solid line with circles in the top two pan- 
els is th e 2-D correlation of t he SDF photometric sample at 
z = 6.6 (|Kashikawa et a'Lll2006l V The insets magnify the region 
4 Mpc/h < r < 20 Mpc/h. 

ticularly on large scales. The impact of cosmic variance is 
relatively independent of the flux sensitivity of the survey. 
Therefore, a larger survey volume than the SDF is necessary 
to mitigate its effect. 

However, Figure [11] suggests that current observations 
in the SDF can distinguish a model with Xi = 0.5 from 
one with Xi ^ 1. Figure [T2l illustrates more explicitly the 
ability to constrain reionization with current and upcoming 
surveys. Note that the y-axis in Figure [12] is linear rather 
than log as in Figure [11] The top two panels in Figure [12] 
show measured from different locations in our box for 
similar survey specifications as in Figure [TT] and with Fc = 
and Fc = 0.25. The dashed curves are 8 randomly selected 
sets of wb for Xi = 0.5, and the solid ones are the same for 
Xi ^ 1. Both groups of curves in the top two panels in Figure 
[12] are measured from independent volumes in the simulation 
box, and each curve is calculated from 58 LAEs. The two 
sets of curves are fairly distinct in the top and middle panels 
in Figure [12] suggesting that Xi < 0.5 can be distinguished 
from Xi ^ 1 using current data. 

The middle panel in Figure [TT] shows the correlation 
function for mock surveys similar to the 2; = 6.6 SXDS that 



will be completed in the next year. We compute this function 
in a 1.4 sq. deg. area (the angular size of the 186 Mpc sim- 
ulation box), with a 130 A narrow band filter, and from 250 
emitters. The SXDS is instead 1 sq. deg., but should have 
a similar number of emitters. Here, the uncertainty in we 
has been reduced substantially compared to the top panel. 
Also see the bottom panel in Figure [T2l With an SXDS-like 
survey, even Xi = 0.7 can, on average, be distinguished from 
Xi ^ 1. 

The bottom panel in Figure [TT] is for a luminosity- 
limited spectroscopic survey in the 1.4 sq. deg. area and with 
Fc = 0, with the same mmin as the middle panel, and with 
the number of spectroscopically confirmed emitters reduced 
by the factor 1 — 0.25 from the mock photometric surveys in 
the middle panel. Since it is more expensive to perform spec- 
troscopy, it is interesting to estimate how much a spectro- 
scopic sample improves the potential to detect reionization. 
The only difference between the top panel and the middle 
panel is that the shot-noise errors have been reduced by the 
factor ^ (1 — Fc). A spectroscopic sample is not significantly 
more sensitive to wb than a photometric sample if Fc — 0.25. 
The benefits of a spectroscopic survey increase for larger 
contamination fractions. In addition, spectroscopic surveys 
open up a third dimension for study. For narrow band LAE 
surveys, the narrow third dimension does not provide much 
additional information, but for broad band Lya surveys (as 
will be conducted with JWST) three dimensional clustering 
statistics will be the way to go. 

Motivated by the appearance of large voids in the LAE 
field, in Appendix [B] we investigate the ability of void statis- 
tics to detect reionization. We find that, provided that 
Fc ~ 0, void statistics can distinguish models that have 
different Xi with comparable significance to the correlation 
function. 

5.3 Subaru Data 

We have seen that a SDF-like data set can distinguish be- 
tween different clustering models for the z — 6.6 LAEs. Let 
us quantify how well the real SDF data can distinguish be- 
tween these different models. We take the z = 6.6 SDF pho- 
tomet ric sa mple of 58 emitter s presented in lTaniguchi et aP 
(2005) and Kashikawa et al.J (j2006l l as well as the exact 
survey specifications, accoun ting for masked region s. We 
then compute we{t) using the lLandv Szalavl (|l993h unbi- 
ased estimator. The connected circles in Figure [T2l represent 
the z — 6.6 cor relation function of emitters presented in 
iKashikawa et al.1 ((2006|). 

We compute the likelihoods of different theoretical mod- 
els for We and Ce given the data (and assuming Gaus- 
sian statistics). These models are constructed in the same 
way as in Section 15.21 We assume a survey that has the 
same angular dimensions as SDF as well as a line-of- 
sight window function W{z) that is a Gaussian in v with 
FWHM 132 A and that is centered at z = 6.6 (such that 
m ey.p{—Ta{'i'Q))W{zm) > mmin iH ordcr to be observed). 
This window function is a fair approximation to the window 
provided by the 9210 A narrow band filter on Suprime-Cam, 
but we find that even if we use a tophat for the window func- 
tion, our conclusions are unchanged, implying that the exact 
functional form of W{z) is unimportant. 

Let us compare the likelihoods of models in which Xi ^ 
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Table 1. = ~2 log£ for different clustering models given the 
clustering in the z = 6.6 SDF photometric sample. These calcula- 
tions assume model (i) for reionization. The model marked "wind" 
assumes = 400 km s"-*^. Note that to achieve the observed 
abundance at z = 6.6 for Xi ^ 1 then, on average, /e = 0.02, 0.1 
and 0.2 for rrimin = 3 x lO^^ M©, 7 x lO^^ M©, and 1 x lO^^ M©, 
respectively. 



mmin (^0) 


Xi 


X [^c — UJ 


X [^c - j) 


3 X 10^0 


1.0 


4.9 


1.1 




0.7 


10.3 


4.9 




0.5 


16.6 


9.9 




0.3 


22.4 


15.1 


7 X lO^o 


1.0 


10.5 


4.9 




0.7 


16.0 


9.6 




0.5 


23.6 


15.9 




0.3 


32.6 


23.9 


1 X lO^i 


1.0 


14.2 


7.8 




0.7 


19.4 


12.2 




0.5 


29.3 


20.7 


3 X 10^° (wind) 


0.7 


8.2 


3.2 




0.5 


11.5 


5.4 




0.3 


13.8 


7.5 



1 to models in which the universe is significantly neutral. 
Table □ summarizes our results. We find that if ?7imin = 7 X 
10^° Mq and Fc = 0.25 (adjusting /e to yield the observed 
tie), a universe with ^ 1 is favored by the data over a 
universe with Xi = 0.5 at 3.3-cr (with Xi = 0.7 at 2.2-cr). For 
a model with mmin = 3 x 10^° Mq and Fc — 0.25 (which 
requires /e = 0.02 when Xi 1), a universe with ^ 1 is 
favored over a universe with Xi = 0.5 at 3.0- a. 

Using reionization model (i), if we marginalize over 
mmin and Fc with the priors mmin > 3 x 10"'^° Mq and 
Fc > 0.25 (with — 0), the maximum likelihood model with 
Xi ^ 1 is favored over the maximum likelihood model with 
Xi < 0.5 at 3.0-a (xi < 0.7 at 1.9-cr/ If we relax the prior to 
mmin > 1 X 10^° Mq (such that at equality /e = 4 x 10~^) 
than Xi < 0.5 is then disfavored at 2-cr. If we used model 
(ii) [model (iii)] to do the calculation discussed in the first 
sentence in this paragraph, a model with ~ 1 is preferred 
over a model with Xi < 0.5 at 2.9-cr [2.8-cr]. Also, for reion- 
ization model (i) and the toy galactic wind model in which 
all LAEs have t'w = 400 km s~^, Xi < 0.5 is disfavored at 
2.1-cr. 

It is useful to note that the clustering data of z = 6.6 
LAEs can even be used to constrain the intrinsic clustering 
of emitters if we assume that x^ ^ 1 at 2; = 6.6. For example, 
the current data favors mmin = 3 X 10^° Mq over TTlmin — 
7 X 10^° Mq at 1.9-cr assuming Fc = 0.25. The current data 
favors Fc = 0.5 over Fc = 0.25 at 1.5-cr for mmin = 3 x 
10^° Mq. 

In this analysis, we have assumed the statistics are 
Gaussian. We use 8 logarithmically-spaced radial bins that 
run between 0.4 and 70 Mpc to compute and Ce- The 
bins at the smallest radii are the least Gaussian because they 
contain the fewest pairs of LAEs. If we discard the first two 
radial bins (r < 1.5 Mpc), the maximum likelihood model 
for Xi ^ 1 can be distinguished by the SDF data from the 
maximum likelihood model for Xi < 0.5 at 2.6-cr rather than 



at 3.0-cr assuming model (i) for reionization and the priors 
mmin > 3 X 10^° Mq, Fc> 0.25, and = 0. 

We ignored the effect of survey incompleteness in the 
above analysis. The SDF survey is complete at ~ 50% level 
above the lumino sity threshold Lmin (and at the ^ 75% level 
above 2.5 Lmin) (|Kashikawa et al.l l2006l V We have investi- 
gated the importance of this effect by randomly discarding 
half of the objects that fall below 1.5 Lmin and one fourth 
of objects between 1.5 Lmin < Lobs,E < 2.5 Lmin, and we 
find that this operation changes the cosmological part of Pe 
by < 25% on all scales and does not significantly alter the 
conclusions in this section. We find that incompleteness is 
degenerate with the parameter /e. 

In this analysis, we have also assumed that the shape 
of the luminosity function was set by the shape of the mass 
function. The z — 5.7 and z — 6.6 luminosity functions 
are poorly constrained at the faint end, providing good 
fit s for faint end power- l aw indexes of /3 = 2, 2.5, and 
3 (|Shimasaku et all l2006l : iKashikawa et all l2006l ). The ef- 
fect of reionization on the LAEs depends on the slope of 
the luminosity function near Lmin because the emitters that 
fall near the detection threshold are the easiest to obscure 
with neutral regions. The luminosity function in our mod- 
els scales as /3 = 2.3 at mmin = 1 x 10^° M©, /3 = 2.6 at 
mmin = 5 X 10^° Mq, and f3 2.8 at mmin = 1 x 10^^ Mq. 
Therefore, the range of mmin that we have considered spans 
much of the relevant parameter space for p. 

The calculations in this section were aimed at under- 
standing the sensitivity to reionization of widefield surveys 
at z = 6.6 that have already been completed or will be 
finished this coming year. These surveys constitute a small 
fraction of the total Subaru observing time in a year. How- 
ever, even with these surveys, constraints can be placed on 
reionization. The sensitivity to reionization of a mission ded- 
icated specifically to this science would be vastly superior. 
We discuss the prospects for a few upcoming LAE surveys 
in Section [6] 



6 FUTURE LAE SURVEYS 

We have shown that if the Universe is largely neutral at 
z < 7, observations like those being done on Subaru should 
be able to study the modulation from the HII regions. What 
is the prospect for upcoming surveys to detect reionization 
at z > 7? To answer this question, we concentrate on three 
instruments: (1) the Dark Ages z Lya Explorer (DAzLE), 
which has started observing on the Very Large Telescope 
(VLT), (2) the configurable multi-slit spectrograph MOS- 
FIRE {the Multi-Object Spectrograph For Infra-Red Explo- 
ration) to be commissioned on the Keck telescope in 2009, 
and (3) NIRSpec on JWST. These represent some of the 
most promising instruments to target high-redshift LAEs. 

Upcoming instruments take different approaches to 
identify high-redshift LAEs. DAzLE uses two 10 A, overlap- 
ping narrow band filters in a 6^ x 6' FoV. Differencing these 
filters makes it possible to subtract out continuum sources, 
leaving just the broad-line sources. The effective volume for 
DAzLE in a single pointing is 1600 Mpc*^, or ^ 1% of the vol- 
ume of the SDF observations. However, DAzLE's flux sensi- 
tivity is below that of other narrow band surveys because of 
its extremely narrow filters. MOSFIRE uses the magnifica- 
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tion bias of cluster lensing caustics to improve the sensitivity 
to high-redshift LAEs by as much as a factor of 50. MOS- 
FIRE is capable of simultaneously taking spectra from 45 
slits, each 7/^3 in length. For an 80 hr observation of 8 clus- 
ters, IStark et al.l (|2007l ) estimate that these surveys could 
cover an effective volume of 80 Mpc^ distributed between 
2; = 7.0 — 8.3- much smaller than the volume of narrow 
band samples. However, whereas DAzLE will be sensitive to 
Lohs,E ~ 10^^ ergs~^ at 2; = 8, MOSFIRE can measure line 
fluxes of lensed galaxies which have Lobs,E ~ 5 x 10^^ erg 
in a 10 hr exposure. 

With the launch of JWST in 2013, the NIRSpec instru- 
ment (capable of simultaneously taking spectra of 100 ob- 
jects in a 3^ X 3^ FoV) will be sensitive to Lya luminosities 
that are over an order of magnitude smaller than terres- 
trial narrow band surveys and comparable to the effective 
sensitivity of a MOSFIRE caustic survey. JWST will use 
the imaging camera NIRCam to select high-redshift galax- 
ies with the Lyman break technique, and follow up with 
NIRSpec to select LAEs. This selection strategy will facili- 
tate cross correlation studies to isolate the bubble-induced 
fluctuations (as discussed in Section [5} and allow 3-D to- 
mography of LAEs. See Figure [3] for a mock survey of LAEs 
that has a JWST-like sensitivity. However, a survey with 
JWST that has the same area as the panels in Figure [3] 
would in fact require over 100 different pointings of both 
NIRCam and NIRSpec to cover this 1000 sq. arc-min. field. 
(The box in the lower left-hand panel of Fig. [3] is the FoV 
for NIRSpec.) The current high-redshift program for JWST 
outlines a deep surve y in a much sm aller area, a tens of sq. 
arc-min. region ( Gardner et al .11200 6). 

In general, it is diflicult to make predictions for the sta- 
tistical significance with which upcoming instruments can 
detect reionization at 2; > 7 owing to all of the uncertain- 
ties in the source properties and their evolution. As we have 
seen, there is substantial uncertainty in the z = 6.6 luminos- 
ity function. The extrapolation of this luminosity function 
to higher redshifts depends sensitively on, for example, the 
typical halo mass of LAEs. If /e = 1 in current LAE surveys 
(such that the observed LAEs are in the most massive ha- 
los), then the luminosity function will evolve more quickly 
with redshift than if /e < 1- Unfortunately, the parameter 
/e is not constrained to even an order of magnitude by the 
z = 6.6 data (Dijkstra et al. 2006; Star k et al. 20071 Other 
factors could also shape the evolution of the Lya luminosity 
function such as the build-up of dust or an evolving abun- 
dance of metal-poor stars within galaxies. 

Despite these challenges, several stu dies have made es- 
timates for tt-e as a function of redshift ("Stark et al.| l2007l : 
iBarton et al.1 l2004l: iThommes M eisenheimer 20o l). For 
DAzL E, iThommes Meisenheimei] (|2005l ) and Star k et all 
(|2007l ) find that roughly 1 LAE should be observed at z = 8 
in a single 10 hr pointing of the camera. To estimate the 
abundance of LAEs that will be observed with MOSFIRE 
on Keck and NIRSpec on JWST is more difficult than for 
DAzLE because these instruments are at least an order of 
magnitude more sensitive, probing a region of the luminosity 
function that has yet to be explored. Fortunately, an esti- 
mate for the observed number of emitters is less important 
for understanding the scientific impact of MOSFIRE and 
JWST because the limiting factor will be cosmic variance 
and not shot noise. 
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Figure 13. Cosmic fluctuations in the measured LAE number 
density as a function of Xi for current and upcoming surveys (the 
variance owing to shot noise has been subtracted). Most curves 
assume one pointing of the relevant imaging camera/spectrograph 
(SDF is one pointing of Suprime-Cam), with the exceptions being 
those for the SXDF and UDS surveys, in which the variance is 
calculated for the survey volume. The DAzLE and Subaru curves 
are computed assuming that mmin = 5 x 10^^ ^0, and those for 
JWST are calculated assuming that mmin = 1 X 10^0 Mq. The 
error bars on the DAzLE and SXDS curves represent the change 
in 6ne/ne if mmin is reduced or increased by a factor of 2. 



Here, we make predictions for the cosmic variance in 
the measured LAE number density for upcoming surveys. 
We do not perform this calculation for a MOSFIRE clus- 
ter caustic surve y owing to the complicated geometry (see 
IStark et alhoO ^f). In the context of our simple model for the 
LAEs, cosmic variance depends on just the parameter mmin- 
When fluctuations in the measured number density of emit- 
ters within a survey volume are on the order of unity about 
the average abundance nE, it will be difficult to understand 
reionization without a larger survey volume. For example, 
these fluctuations will prevent a precise measurement of the 
luminosity function. 

The curves in Figure [13] represent the standard devi- 
ation in the number density of LAEs (with the shot noise 
component subtracted out) measured with current and up- 
coming surveys (and assuming no foreground contamina- 
tion) . Note that neutral gas increases the standard deviation 
signiflcantly. The DAzLE and Subaru curves are computed 
assuming that mmin = 5 x 10^° Mq (which is a conserva- 
tive choice and requires /e < 0.1), and the ones for JWST 
are calculated assuming mmin = 1 x 10^° Mq. All of these 
curves are determined by first convolving the LAE maps 
with a window function that has the same volume specifica- 
tions as these surveys/instruments and, then, by computing 
the variance of the windowed field. For JWST, which will 
take spectra over a broad range of wavelengths with NIR- 
Spec, we have assumed that the line-of-sight width used to 
estimate nE{z) is 100 Mpc {Az ^ 0.3). The choice of 100 
Mpc is motivated by the short duration of reionization in 
simulations of reionization (0.5 < < 1 spans Az ~ 1) and 
the desire to study LAEs for a few different ionized fractions 
during reionization. 

The curves in Figure [13] are computed from snapshots 
between z = 7 and z = 8. The value of Sue /tie depends 
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modestly on the redshift of the survey. For example, for 
the SXDS, for mmin = 5 x 10^° Mq and Xi ^ 1 then 
SuE/riE — 0.13 dX z — 7 whereas 5n^/n^ — 0.17 dX z — 80 
In addition, if we vary the value of mmin, the results do 
not change significantly. The error bars on the DAzLE and 
Subaru-SXDF curves in Figure [13] represent the change in 
dne/rie if the fiducial mmin is reduced or increased by a fac- 
tor of two. 

The volume surveyed by DAzLE is much smaller than 
the volume of a 20 Mpc bubble. For some pointings, the 
surveyed volume may consist entirely of a neutral region 
(and zero emitters will be observed), and, for other point- 
ings, it may fall entirely within an HII region. Figure [T3l 
shows that in order to constrain the cosmic fluctuations in 
the luminosity function to 50%, DAzLE needs Np ^ 10 
non-contiguous pointings when Xi ^ 1 and significantly 
more pointings when the Universe is largely neutral. (Note 
that (SnE)Np ~ 5n^/ y/N^ for non-contiguous pointings.) Of 
course, this observing strategy is optimized for reducing the 
cosmic variance. Shot noise may be as important for DAzLE. 

With JWST, one might hope to be able to place a tight 
constraint on the LAE luminosity function. However, for 
a 20% constraint, JWST needs Np ^ 10 non-contiguous 
pointings when ~ 1 (and Np ^ 50 when Xi ~ 0.5). 

Much of the high-redshift data for JWST will be gath- 
ered in the Ultra-Deep Survey (UDS), which will be a con- 
tiguo us field spanning tens of sq. arc-min. ([Gardner et al.l 
I2OO6I I. Here, we generously assume a square field of 80 sq. 
arc-min 9 tiles observed with both NIRCam and NIR- 
Spec). Figure [13] shows that this observation is far from opti- 
mal for constraining the luminosity function with 5n^/n^ > 
0.3. A different observing strategy is necessary for JWST to 
derive tight constraints on reionization with LAEs. 

The largest volume LAE survey currently being con- 
ducted at 2; > 8 is the MOIRCS Deep Survey on Subaru, 
which targets LAEs at z — 8.8. While not as sensitive as 
DAzLE, its volume of 7' x A' x 10 A is an order of magni- 
tude larger than that of DAzLE (jOuchi et al.ll2007l ). If we 
had plotted a curve of 5n^/n^ for MOIRCS in Figure [13] 
this curve would fall between the curves for the Subaru-SDF 
and for DAzLE. 

In this section, we have seen that cosmic variance will 
be a significant concern for upcoming LAE surveys at 2; > 7. 
Degree- scale surveys make the cosmic variance manageable 
for this science. 



7 CONCLUSIONS 

LAE surveys are probing increasingly higher redshifts. Not 
only will these surveys inform us about high-redshift galax- 
ies, but they have the potential to be the first observations 
to unambiguously detect patchy reionization. The tens of 
Mpc HII regions during reionization modulate the observed 
distribution of LAEs and boost their observed clustering. 
We have shown that this effect on the angular correlation 
function (Section[5]) or the void probability distribution (Ap- 
pendix [BJ can be well in excess of the intrinsic clustering of 

Note that the SXDS curve in Figure [T3l assumes a square sur- 
vey in a 1 sq. deg. FoV. SXDS is instead a cross-shaped survey 
with the same volume. 



halos in the concordance cosmology. This enhanced cluster- 
ing depends most strongly on x^, but also somewhat on the 
morphology of the HII regions during reionization. Observ- 
ing enhanced clustering would confirm the prediction that 
the HII regions during reionization are large. 

Even the current 2; = 6.6 Subaru LAE survey, which 
has a photometric sample of only 58 LAEs in a 0.25 deg^ 
field, can place constraints on the reionization process. We 
find that the angular correlation function of the SDF pho- 
tometric sample of z — 6.6 LAEs favors an ionized uni- 
verse over a universe with Xi < 0.5 at a 3-cr confidence 
level (2-cr if all emitters have strong galactic winds). This 
constraint is both competitive with oth er constraints on Xj 
obtained from GRBs and from LAEs (jTotani et al.1 12OO6I : 
iMalhotra Rhoadsll2006bV and it rules out the picture of 
the z = 6 quasars expanding into a neutral IGM. In addi- 
tion, this is the first constraint on the ionized fraction that 
consistently accounts for patchy reionization - the favored 
model for how the universe is ionized. Observations in the 
next year in the SDXS will increase the z = 6.6 sample by 
a factor of 5 and place even stronger constraints on Xi. 

Detecting reionization through LAE clustering also of- 
fers a simple consistency check for whether the observed 
correlations owe to reionization: observe the same field with 
a second selection technique in addition to selecting objects 
by their Lya emission. If the Lya-selected galaxies show en- 
hanced clustering relative to the galaxies selected with the 
other technique, the evidence for reionization is strength- 
ened. Similarly, a comparison of the phases of Fourier modes 
between the two galaxy samples can be used to detect reion- 
ization. 

While measuring enhanced clustering is the most fool- 
proof method to detect reionization with LAEs, the presence 
of neutral regions in the IGM also influences the properties 
of the Lya line and the luminosity function. Observing evo- 
lution that is consistent with reionization occurring in all of 
these different statistics would strengthen the argument for 
reionization. We have made predictions for the scale of these 
effects. We show in Appendix [A] that the effect of reioniza- 
tion on the average line profile is less significant than most 
previous studies have found. This result owes to the much 
larger HII regions that arise from properly treating source 
clustering and because emitters that are not significantly ob- 
scured by the neutral regions are preferentially observed. Be- 
cause of these effects, combined with the uncertainty in the 
astrophysical processes that determine the Lya Une shape, 
the consequences of reionization will be difficult to isolate in 
the line profile. 

The impact of reionization on the luminosity function is 
a more promising diagnostic than the line profile. Its effect 
on the luminosity function is similar among the three reion- 
ization models we have considered, if we compare at fixed Xi. 
In all models and at all Xi, the LAE luminosity function is 
suppressed by a factor that is fairly constant as a function of 
luminosity. Our predictions for the impact of reionization on 
the luminosity function are fairly robust to our assumptions 
concerning the bias and intrinsic luminosity of emitters, but 
it may be difficult to distinguish evolution in the luminosity 
function that owes to reionization from evolution owing to 
changing intrinsic properties of the emitters. 

It is likely that most of reionization occurs at higher 
redshifts than probed by current samples. Several LAE sur- 
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veys will target z > 7 LAEs in the coming years. Unfortu- 
nately, these upcoming studies are not optimal for measuring 
the LAE luminosity function or for detecting reionization- 
induced clustering, having fields of view that are too small 
to measure this effect precisely. Even the multi-billion dollar 
satellite JWST has this design flaw. JWST's enhanced sen- 
sitivity over that of current telescopes still makes it useful 
for this science. However, a multi-billion dollar space mis- 
sion is unnecessary to understand reionization with LAEs. 
Widefield observations like those on Subaru, but that target 
z > 7, will be able to put constraints on the reionization 
process. 
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APPENDIX A: LINE PROPERTIES 

Several studies have attempted to understand the impact 
that reionization has on the profil e of Lyg l i nes using sim- 
ple analytic niodels (Miralda-Escude [199^ : ISantosI l2004l: 
Haiman k Cenll2005l : iFurlanetto et al.ll2004 j : IPiikstra etaP 
2OO7I I. [n this section, we investigate the effect of reioniza- 
tion on the line profile using simulations. We show that even 
when a significant fraction of the universe is neutral, damp- 
ing wing absorption owing to neutral patches in the IGM 
tends to have a small effect on the average line shape. We 
use Method 1 described in Section 12.21 to calculate the in- 
trinsic line profile. 

For simplicity. Method 1 sets the intrinsic width of the 
line to be Aiy = VaVvwIc^ where v^w is the circular velocity 
of the host halo at the virial radius. This Az/ is the width 
one expects if the typical column density of HI inside the 
galaxy, A^hi, is low enough such that ionizing photons can 



escape from the galaxy (A^hi ^ 10 cm ) and, for the mas- 
sive galaxies of interest, ev en if the colunin den sity is as high 
as AThi = 10^^"^°cm-2 ( Diikstra et alll2007l ). This Az/ is 
also consistent with the emergent line profile of a 2-phase 
ISM in which dense clouds of HI (with dust inside) scatter 



Lya photons. A 2-phase structure may be necessar y for the 
transmissi on of Lyg photons in the presence of dust (Ne ufeldl 
[1991; Hansen h Oh 2006). We do not attempt to parame- 
terize all of the complicated processes that determine the 
width and shape of the intrinsic line, such as the geometry 
of the galaxy and its bulk motions. Such processes can add 
additional features to the line profile and make it even more 
difficult to disentangle the effect of reionization. 

In our calculations, the observed line profiles differ from 
one another primarily because of two effects: (1) The extent 
along the line of sight from the emitter to the HII bubble 
edge. The extent of the HII region determines the strength of 
the damping wing absorption. This effect is the most impor- 
tant for our study. (2) The amount of infall of gas around the 
emitter. This dictates how much of the red side of the line 
is absorbed as some photons on the red side must redshift 
through resonance because of the Doppler effect, typically 
resulting in absorption (Outflows - which we incorporate 
with a simple prescription - have the opposite effect, allow- 
ing more of the blue side of the line to be transmitted.) The 
effect of infall is only crudely modeled here because these 
computations are performed on a 0.36 Mpc grid, whereas 
the virial radius of a 5 x 10^° M© halo is 0.09 Mpc at 2; = 7. 
This results in an under-prediction of the effect of infallF^ 

The quality of the spectrum from a single 2; > 6 LAE 
for current ground-based telescopes is never good enough 
to constrain in detail the shape of a single Lya line. Pre- 
vious LAE surveys have averaged the line profile from all 
of their emitters to generate a higher quality effective line 
profile. Figure [ATI shows the average line profiles calculated 
from snapshots with Xi — 0.3 (io'p panel), Xi = 0.5 {middle 
panel), and Xi = 0.7 {bottom panel) for z = 8.2, z = 7.7 
and 7.3, respectively. We have normalized the integral over 
these curves to unity. The thick dashed curves are the aver- 
age profile for LAEs with 0.5 Lint, e(2 x lO^^M©) < Lobs,E < 
0.5 Lint, e(4 X 10^^ Mq) for the quoted Xi, whereas the thick 
solid curves are for the same luminosity range but with 
Xi — Because of the exponential falloff in the halo mass 
function at large m, if we eliminated the upper bound on 
Lobs,E and included more luminous LAEs in the average, 
the line profiles would not be significantly affected. 

Note that the differences between the dashed and solid 
curves in Figure [ATI are small, particularly for the Xi = 0.5 
and Xi = 0.7 cases. The reason why reionization does not 
have more of an effect on the line profiles is because the av- 
erage bubble size in which observed LAEs lie is larger than 



■"^^ The blue side of the line is almost always entirely absorbed 
in our calculations. While most of the blue side of the line will 
typically be absorbed in reality, our calculations overestimate this 
effect because the simulations do not resolve the scale on which 
pressure smooths out density fluctuations. This lack of resolution 
artificially suppresses the transmission. The resolution is sufficient 
to capture the most important effect for this paper, the damping 
wing absorption. 

Analytic studies typically assume that the infall is maximal 
at the viral radius and equal to t he circular velocity of the halo 
(|Santosll2004l : iDiikstra et al.ll2007l) . Galactic winds can collision- 
ally ionize the hydrogen within up to ~ 0. 2 pMpc of a galaxy, 
eliminating the importance of nearby infall (Santos 2004). 

The 0.5 that appears in these bounds is there to approximate 
the effect of resonant absorption, which typically absorbs all fre- 
quencies blue- ward of uq. 
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Figure Al. Observed line profile computed from simulation out- 
puts in which xi = 0.3 (top panel), Xi = 0.5 (middle panel), and 
Xi = 0.7 (bottom panel) for z = 8.2, z = 7.7, and z = 7.3, respec- 
tively. All curves have been normalized such that their integral 
is unity to facilitate comparison of the line profile (see Fig. IA2I 
for the relative normalization). The thick dashed curves are the 
average profile of LAEs with 0.5 Lint e(2 x lO^^ Mq) < L^bs E < 
0.5 Lint, e(4 X 10"*^° ^o) for the quoted Xi whereas the thick solid 
curves are the same for the post-reionization, ionized-IGM case. 
The thin dashed and thin solid curves are the same as the thick 
curves but for the LAEs with 0.5 Lint, e(1 x 10^^ Mq) < Lobs,E < 
0.5 Lint, e(2 X 10^^ Mq). In the top panel, the dotted and dot- 
dashed curves that are offset to the right are respectively the the 
Xi 1 and the Xi = 0.3 average spectrum for a toy wind model 
discussed in the text. These calculations are performed using the 
94 Mpc simulation of Model (i). 

a pMpc for the x^ shown in Figure [AT] Previous calculations 
assumed that each LAE was responsible for ionizing its own 
bubble, resulting in i^^ < 1 pMpc (eqn. [2]) and, therefore, 
a larger effect of the damping win ^ absorption on the ob- 
served line profile (e.g., Santos|200j|). One might wonder why 
reionization has any affect on the LAE luminosity function 
and correlation function (as seen in the body of this paper) 
since the average line profile is not strongly affected. Even 
though the line profile is not changed, the average transmis- 
sion is decreased by reionization, and many emitters are not 
observed because of reduced transmission. 

Galactic winds/outflows, which are found to be very 
prevalent in lower redshift LAEs, allow more of the blue 
side of the line to be transmitted. Here we adopt a simple 
toy model for galactic winds which is meant to exaggerate 
their potential effect. We make the same assumptions for the 
intrinsic line profile as before, but we redshift the entire line 
by the velocity of the wind, which is taken to be 400 km s~^ . 
This value for is motived by the average velocity offset 
of th e Lya line in strong emitters at 2; ^ 3 (Shaple v et al.l 
I2OO3I ). Since galaxies are less massive at z > 6 than at 2; ~ 3, 
the galaxies probably cannot drive such powerful winds and 
= 400 km s~^ can be thought of as an upper bound on 
v^. Furthermore, it could be the case that only a fraction of 
the Lya photons are scattered by the wind. The difference 
between dotted and dot-dashed curves in Figure lAll illus- 
trate the effect of reionization in this wind model. Winds 
reduce the importance of damping wing absorption. 

Current narrow band surveys most likely have spectro- 
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Figure A2. Ratio of the average observed luminosity during 
reionization to the average observed luminosity for an ionized uni- 
verse as a function of halo mass. The dot-dashed, dotted, and solid 
thick curves are computed from simulation snapshots of Model (i) 
which have Xi = 0.3 (z = 8.2), x^ = 0.5 (z = 7.7), and Xi = 0.7 
(z = 7.3), respectively, and assuming a survey that is sensitive 
to Lobs,E > 0.5 Lint, e(1 X 10-*^° ^o)- These curves are computed 
with Method 1 in Section [2] The thin curves are the average of 
exp(— ra(z^o)) (as is computed for Method 2 in Section O, using 
the emitters that are observed in the calculation with Method 1. 
The similarity between the thin and thick curves justifies the use 
of Method 2. 

scopically confirmed LAEs that have m ~ 10^^ Mq rather 
than m ^ 2 X 10^° Mq, as considered above. (Only the 
brightest emitters in a narrow band survey can be spec- 
troscopically confirmed.) If we take instead emitters within 
the range 0.5Lint,E(l x 10^^ Mq) < Lint,E < 0.5Lint,E(2 x 
10^^ Mo) (thin curves in Fig. \A1^ , the differences between 
the line profiles for x^ ~ 1 and Xi ~ 0.5 are even smaller than 
in the case previously examined, owing to these more biased 
sources sitting preferentially in the largest HII regions. 

In our computation, the FWHM of the observed line 
profiles for LAEs that have m ^ 2 x 10^° Mq is ^ 
4 A (thick curves in Fig. \A1^ and for those that have 
m ^ 10^^ Mo is ^ 8 A (thin curves). The FWHM of 
the average line of the spectroscopically confirmed z = 
6.6 emitters is ^ 10 A, slightly larger than the FWHM 
in our calculatio ns ( Kashikawa et al. 2006). Interestingly, 
iKashikawa et al.1 (|2006l ) finds a weak ant i- correlation be- 
tween the FWHM of the line and Lobs,E. In our simple 
model for the LAEs, in which Au ^ m^^^, there should be 
a correlation. iHaiman Ce^ (|2005l l suggest that an anti- 
correlation may b e a signature of r eioniz ation. However, to 
reach this result, iHaiman CenI (|2005h assume that the 
intrinsic widths of the Lya lines are the same for emit- 
ters of all Lint,E and that the HII region around an emit- 
ter is cre ated just by this emitt er. The ant i- correlation re- 
sult from IHaiman CenI (|2005l ) depends on these dubious 
assumptions. There is probably too much uncertain astro- 
physics to understand reionization through the correlation of 
the FWHM with Lobs,E- The weak ant i- correlation between 
FWHM with Lobs,E that is observed (with low statistical 
significance) may indicates that there is more dispersion in 
the mapping between m and Lobs,E than we have assumed 
in this work or that winds are important in some z = 6.6 
emitters. As discussed in previous sections, increasing the 
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dispersion in the mapping between m and Lobs,E does not 
significantly affect our conclusions pertaining to the LAE 
luminosity function and to LAE clustering. 

The emitters in a survey whose line profiles are most 
strongly affected by reionization may constitute a small sub- 
sample of the LAEs in a survey. Therefore, it might be more 
fruitful to look at individual lines rather than the average 
line profile to detect reionization. However, even if a survey 
has sufficient signal to noise to study the Lya line from a 
single emitter, it will be challenging to isolate the impact of 
damping wing absorption from other effects. 

Thus far, we have considered the effect of a single model 
for reionization, model (i), on the line profile. In model (iii), 
since the bubbles are smaller, the average line profile at fixed 
Xi is slightly more altered by reionization than in model (i). 
The opposite is true in model (ii), in which the bubbles are 
larger. However, the differences in the line profiles between 
the models at fixed Xi are comparable to or smaller than the 
differences between the line profiles at Xi — 0.3, 0.5, and 0.7 
in model (i). 

Figure IATI ignored the effect of reionization on the trans- 
mission of the Lya line and just illustrated its effect on the 
line profile. The thick curves in Figure [A2l represent the ra- 
tio of the average of Lobs,E at Xi to the same but at ^ 1. 
Only emitters with Lobs,E > 0.5 Lint, e(1 x 10^° ^o) are 
included in these averages. This ratio yields approximately 
the factor exp(— Tc(z^o)) - the strength of damping wing ab- 
sorption. The thick dot-dashed curve is for Xi = 0.3, the 
thick dotted curve is for Xi = 0.5, and the thick solid curve 
is for Xi = 0.7. As Xi decreases, the average Lobs,E also 
decreases owing to smaller bubbles and increased damping 
wing absorption. At low m, the decrease in Lobs,E with x^ is 
smaller than at large m. This result might may appear sur- 
prising since more massive halos sit in larger bubbles, which 
would result in the opposite tendency. However, this trend 
is simply a selection effect owing to the luminosity threshold 
of this mock survey. In addition, this selection effect illus- 
trates why the line profiles are not affected as substantially 
as one might expect: Owing to the nature of the luminosity 
function, most emitters sit near the luminosity threshold. 
In order for these emitters to be observed, they must sit in 
very large HII regions such that the effect of damping wing 
absorption is small. 



APPENDIX B: VOID STATISTICS 

The power spectrum may not be the optimal statistic to 
pick out the modulation of LAEs owing to the bubbles. One 
distinctive feature of the LAE field during reionization is 
that there are large voids owing to the large neutral regions 
(see Figs. [2] and [3]). This feature motivates using the void 
probability distribution function P (R, 0) to measure the e f- 
fect of reionization, as was done in iKashikawa et all (|2006h . 
Void statistics have a long history of being used t o interpret 
galaxy surveys (e.g.. White 19791 and ICroton ~ al. 2004). 

P{R, 0) is the probability that a circle of radius R 
around a point in the survey field does not contain any 
galaxies. For a Poisson distribution of galaxies, P(i?, 0) = 
exp(— 7ri?^ ^e), where Ee is the surface density of LAEs. A 
generalization of this statistic is the probability that a region 
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Figure Bl. 2-D void probability distribution function P(R, 0) at 
z = 7.5 computed from a simulation of model (i). Each curve is 
calculated by averaging five independent surveys of 200 emitters 
in a volume of 186 Mpc x 186 Mpc x 130A. The error bars are 
the standard deviation on the curves, calculated from five mock 
surveys in independent volumes within the simulation box. The 
thin solid curve is P{R, 0) for a purely Poisson distribution. The 
thick solid and dot-dashed curves are P(R, 0) for Xi ^ 1 with 
LAEs in halos with m > 5 x 10^° M© (/e = 0.15) and with 
m > 1 X 10^^ Mq (/e = 0.6), respectively. The long dashed [short 
dashed] curves are the ionized case for Xi = 0.5 and m > 3 x 
10^0 Mq (/e = 0.15) [for Xi = 0.5 and m > 5 x lO^^ Mq (/e = 
0.5)]. 
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Figure B2. The same as IBll but the curves are calculated at 
z = 6.6 and rrimin = 5 x 10^^ ^Q- This figure compares the 
P(R,0) in reionization model (i) for different Xi. 



contains N galaxies, P{R,N), and, for a Poisson distribu- 
tion, P{R,N) = (ttR^I^e)^ exp(-7ri?2 Se)/A^!. 

Figure [BT] displays P(i?, 0) computed from the 186 Mpc 
simulation of model (i). Each curve is calculated by averag- 
ing five independent surveys of 200 emitters in a volume of 
186 Mpc x 186 Mpc x 130A (a little larger than the volume 
of the SXDS), yielding ue = 1.5 x 10"^ Mpc-^ The thin 
solid curve is P{R, 0) for a Poisson distribution. The thick 
solid and dot-dashed curves are P(i?, 0) of an ionized uni- 
verse with LAEs in halos with m > 5 x 10^*^ Mq and with 
m > 1 X 10^^ respectively. The dashed (dotted curves) 
are the ionized case with Xi — 0.5 and with m > 3 x 10^° Mq 
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(m > 5 X 10^° Mq). We adjust the duty cycle to fix the 
number of emitters, requiring /e = 0.5 for the dotted curve, 
/e = 0.15 for the dashed curve, /e = 0.6 for the dot-dashed, 
and /e = 0.15 for the thick solid. 

The error bars in Figure IBll are the standard devia- 
tion from a sample of five independent mock surveys with 
the specifications given above. Note that the error bars in 
different bins are correlated. From the errors, we see that 
we can distinguish the models in this plot at > 1-a confi- 
dence level, and, in particular, we can distinguish the case 
with reionization from those without reionization. This is 
promising because the dot-dashed curve assumes that the 
LAEs are in the rarest, most clustered halos, and this curve 
is below those with Xi = 0.5. Clustering of LAEs owing to 
reionization generates large-scale voids to a much larger de- 
gree than the intrinsic clustering of halos in ACDM. 

Figure [B2l shows P{R, 0) for different Xi calculated from 
reionization model (i). These curves are calculated assuming 
200 LAE, mmin = 5 X 10^*^ ^o, and for a survey compara- 
ble in volume to the SXDS. As with the correlation function, 
P(i?, 0) can distinguish the four different Xi in this figure. 
Qualitatively, the significance with which P(i?, 0) enables 
one to distinguish between these different Xi appears to be 
comparable to the significance the correlation function af- 
fords. 

Thus far, we have assumed no contamination from fore- 
ground galaxies, which is certainly not the case for pho- 
tometric LAE surveys. Unfortunately, a luminosity- limited, 
widefield spectroscopic LAE survey at ^ > 6 has not 
been conducted. If there is significant contamination in 
the survey, then the measured value of P{R,0) will be 
significantly biased, suppressed (on average) by the factor 
exp(— 7ri?^ Scont) where Scont is the density of contaminat- 
ing galaxies (assuming that the contaminants are uncorre- 
lated). If Scont/S is appreciable, t his suppress i on is signif- 
icant, especially at large R. Kashikawa et al.l (|2006h mea- 
sured P{R, 0) from the SDF photometric sample of LAE 
at z = 6.6. They found that the measured P{R,0) is con- 
sistent with a random sample. This conclusion may owe to 
contamination biasing their measurement. 

We can attempt to alleviate this issue of foreground 
galaxy contamination biasing P(i?, 0) by instead computing 
the statistics P{R, 1), P(i?, 2), etc., which are progressively 
less biased by this effect. In this vein, we can do the exact 
opposite and use peaks rather than voids to probe reioniza- 
tion. As long as these peaks are sufficiently rare, they will 
exist only in a distribution of emitters that is clustered. In 
addition, peaks in the emitter distribution are a distinctive 
feature of Lya maps during reionization (see Figs. [2] and [3|). 

We have investigated the statistic 

oo 

PAR,y)= J2 PiR,n), (Bl) 

n=floor(7 ttR'^ S) 

where 7 is a constant that we have varied between 3 and 
10 and floor(x) is a function that returns the largest inte- 
ger smaller than x. We find that Pp{R,j) can distinguish 
the models considered in Figure [BTI with comparable signif- 
icance to P(i?, 0). Like P(i?, 0), a survey needs an estimate 
for its contamination fraction to be able to predict Pp(R,j). 
This can be estimated by performing spectroscopic follow- 
up on a portion of the survey. Unlike P(i?, 0), the bias in 
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Figure CI. Effect of the duty cycle of ionizing sources on the 
morphology of reionization. The white regions are ionized and 
the black are neutral. The left panels are from a simulation using 
model (i), and the right panels are from an identical simulation 
except that only 10% of the sources are active at any time, and 
these sources have 10 times the luminosity compared to the same 
sources in model (i). The active sources are chosen randomly, and 
this randomization is performed every 20 million years. The top 
panels have Xi = 0.3, and the bottom panels have Xi = 0.7. All 
panels are from a slice through the 186 Mpc box. 

the estimation of Pp(i?,7) incurred by foreground contami- 
nation is not an exponential factor, and, therefore, a survey 
can more accurately correct for this bias with an estimate 
of Pc. 

More work needs to be done to quantify how much inde- 
pendent information derives from void (and peak) statistics 
compared to the correlation function and to quantify the 
benefits of these statistics over the correlation function. 



APPENDIX C: DUTY CYCLE 

In iMcQuinn et al.l (|2006l l , the effect of the morphology of 
reionization on the source properties, on the redshift of 
reio nization, and on theminihalos were investigated in de- 
tail. iMcQuinn et~al] (|2006h found that the structure of the 
HII regions during reionization was robust to most consid- 
ered effects, with the most important dependence being the 
value of the ionized fraction. However, the effect of a duty cy- 
cle for the ionizing sources - the fraction of time the sources 
are emitting ionizing photons - on the structure of the HII 
regions was omitted in their analysis. It is probable that 
star formation at high redshifts is sporadic. At any given 
time, a small fraction of galaxies may contain the high mass 
stars that produce ionizing photons. These active galaxies 
are probably also the LAEs. 

It is interesting to understand the effect of a small ion- 
izing photon duty cycle /ion <C 1 on the structure of the HII 
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regions. None of the simulations used in the body of this pa- 
per included this effect. Figure [CTI compares the 130 Mpc/h 
simulation using model (i) in which /ion = 1 {left panels) 
to a simulation where /ion = 0.1 and in which the sources 
that are active have a luminosity that is boosted by a fac- 
tor of 10 over the sources in model (i) (right panels). To 
achieve /ion = 0.1, 10% of the sources are randomly selected 
to be active, and this randomization is performed every 20 
million years. Twenty- million years was chosen to roughly 
match the lifetime of massive stars. One can see in Figure 
ICll that the large HII regions are very similar between these 
two cases. This invariance owes to thousands of emitters in 
each large bubble contributing to its ionization, such that 
the total number of photons produced inside the bubble, 
and, therefore, the bubble size is largely unaffected by /ion- 
A similar conclusion about the effect of /ion was reached in 
the analytic study of iFurlanetto et al.l (120051 ). 

The duty cycle does influence the smallest HII regions 
- bubbles small enough such that shot noise in the number 
of galaxies is important. The criterion for /ion to change 
the LAE statistics from models in which all galaxies are 
active is that the fluctuations owing to shot noise in the 
total number of source galaxies must be comparable to the 
cosmological fluctuations in the number of these galaxies 
within bubbles of size > 1 pMpc. At z ^ 7 there are ^ 
10"^ /ion galaxies emitting ionizing photons in an HII region 
of size 1 pMpc in the source models discussed in this paper, 
whereas cosmological fluctuations in the source abundance 
are of the order of unity at this scale. Therefore, it would 
require either /ion 10~^ for the value of /ion to affect the 
statistics of the LAEs. 



cross correlated the ionization field generated from a simula- 
tion that has source grouping to one without source group- 
ing, but with the same sources and at fixed Xi. We find no 
appreciable diflFerence between the two ionization fields, im- 
plying that source grouping does not affect our conclusions. 
Source grouping speeds up the radiative transfer algorithm 
considerably when Xi > 0.5. 

More aggressive grouping does lead to less accuracy in 
determining the local photo-ionization rate around sources, 
which can infiuence the transmission properties of the Lya 
line. However, we only use this more aggressive grouping 
algorithm in combination with Method 2 in Section [2] to 
calculate the fine fiux, and Method 2 does not depend on 
the details of the photo-ionization state of the gas. 



APPENDIX D: SOURCE GROUPING 

The lSokasian et"aL 1 (|200lL I2OO3I I and lMcQuinn et al.1 (|2006l l 

radiative transfer code has been improved for this paper to 
group sources more efficiently, speeding up the computation. 
In previous versions of the code, only sources that fall within 
the same grid cell were grouped as a single source, often 
resulting in tens of millions of sources that must be processed 
in a time step. This algorithm did not take advantage of the 
fact that the HII regions become much larger than the size of 
a grid cell, motivating a more aggressive grouping algorithm 
without a significant loss in accuracy. 

Our algorithm for grouping sources is as follows: 

(i) Smooth ionization field Xi at scale R with a top hat 
filter, yielding the field Xi,R. 

(ii) Loop over sources in random order. 

If a source is in a cell in which Xi^R > check that 
Xi,r > C is also satisfied at all r < i?. If yes, group sources 
within a sphere of radius r]R into a single source, and place 
the new source at the center of luminosity of the grouped 
sources. 

(iii) Repeat previous steps for smaller R, but only using 
the sources which have not been grouped. Stop when r/ R is 
less than the width of a grid cell. 

(iv) Take the new grouped source field and repeat all pre- 
vious steps once again. 

In this paper, we start smoothing with = 50 Mpc, and set 
C = 0.9 and 77 = 0.1. As a test of this algorithm, we have 
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